Logo Search packages:      
Sourcecode: t-coffee version File versions  Download package

tree_util.c

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <stdarg.h>
#include <strings.h>

#include "io_lib_header.h"
#include "util_lib_header.h"
#include "dp_lib_header.h"
#include "define_header.h"

#define TOPOLOGY 1
#define WEIGHTED 2
#define LENGTH   3
#define RECODE   4

int distance_tree;
int rooted_tree;
int tot_nseq;
/*********************************************************************/
/*                                                                   */
/*                                   dpa_tree_manipulation           */
/*                                                                   */
/*                                                                   */
/*********************************************************************/
static NT_node code_dpa_tree ( NT_node T, int **D);
NT_node collapse_sub_tree ( NT_node T,int nseq, int *list, char *new_name);

NT_node tree2dpa_tree (NT_node T, Alignment *A)
{
  /*This Function sets the branches with Length values used by DP*/
  /*The tree must be rooted*/
  Sequence *S;
  int **D;
  
  
  S=aln2seq (A);
   
  T=recode_tree     (T, S);
  T=reset_dist_tree (T, -1);
  D=get_sim_aln_array (A, "idmat");
 
  
  T=code_dpa_tree (T, D);
  return T;
}

NT_node code_dpa_tree ( NT_node T, int **D)
{
  if ( !T) return T;
  else if ( T->leaf==1)
    {
      T->dist=100;
      return T;
    }
  else
    {
      int nl, *ll;
      int nr, *lr;
      int a, b, min=100;
      float tot, n=0;
      
      nl=(T->left)->nseq;ll=(T->left)->lseq;
      nr=(T->right)->nseq;lr=(T->right)->lseq;

      for (tot=0,n=0, a=0; a< nl; a++)
      for ( b=0; b< nr; b++, n++)
        {
          tot+=D[ll[a]][lr[b]];
          min=MIN(min,(D[ll[a]][lr[b]]));
        }
      /*      T->dist=(mode==AVERAGE)?(tot/n):min:;*/
      T->dist=(n>0)?tot/n:0;
      T->dist=min;
      code_dpa_tree ( T->right, D);
      code_dpa_tree ( T->left, D);
      return T;
    }
}
static int group_number;
char *tree2Ngroup (Alignment *A, NT_node T, int max_n, char *fname)
{
  double top, bot, mid, pmid;
  Sequence *S;
  int n;
  
  

  if (!T)
    {
      
      T=tree_compute(A, 0, NULL);
      T=tree2dpa_tree (T,A);
    }
  S=tree2seq(T, NULL);
  
  if ( max_n>S->nseq)max_n=S->nseq;
  
  

  top=100; bot=0;
  pmid=0; mid=50;
  n=tree2group_file(T, S,0, (int)mid,fname);
  mid=dichotomy((double)n, (double)max_n,(pmid=mid), &bot, &top);
  while (n!=max_n && pmid!=mid)
    {
      n=tree2group_file(T, S,0, (int)mid, fname);
      mid=dichotomy((double)n, (double)max_n,(pmid=mid), &bot, &top);
    }
  fprintf ( stderr, "\n#TrimTC: Split in %d Groups at a minimum of %d%% ID\n",n, (int)mid);
  return fname;
}  
static int group_number;
int tree2group_file ( NT_node T,Sequence *S, int maxnseq, int minsim, char *name)
  {
    FILE *fp;
    int n;
    
    fp=vfopen (name, "w");
    vfclose (tree2group (T, S,maxnseq,minsim, "tree2ngroup",fp));
   
    return count_n_line_in_file(name);
  }
    
      
FILE * tree2group ( NT_node T,Sequence *S, int maxnseq, int minsim,char *name, FILE *fp)
{
  if ( !T)return fp;
  else
    {
      int m,d;

      m=(maxnseq==0)?S->nseq:maxnseq;
      d=minsim;


      
      if ( T->nseq<=m && T->dist>=d)
      {
        int a;
        fprintf ( fp, ">%s_%d ", (name)?name:"", ++group_number);
        for ( a=0; a< T->nseq; a++)
          fprintf ( fp, "%s ", S->name[T->lseq[a]]);
        fprintf (fp, "\n");
        if (!T->parent)group_number=0;
        return fp;
      }
      else
      {
        fp=tree2group (T->right, S, maxnseq, minsim, name,fp);
        fp=tree2group (T->left, S, maxnseq, minsim, name,fp);
        if (!T->parent)group_number=0;
        return fp;
      }
      
    }
}


NT_node  tree2collapsed_tree (NT_node T, int n, char **string)
{
  char ***list;
  Sequence *A;
  int a, *nlist;

  
  A=tree2seq(T, NULL);
  T=recode_tree(T, A);
  list=vcalloc (A->nseq, sizeof (char***));
  nlist=vcalloc (A->nseq, sizeof (int));
  if ( n==0)return T;
  else if (n>1)
    {
      int l;
      char *buf;
      
      for (l=0,a=0; a< n; a++)l+=strlen (string[a]);
      buf=vcalloc ( 2*n+l+1, sizeof (char));
      for (a=0; a< n; a++){buf=strcat (buf,string[a]), buf=strcat ( buf, " ");}     
      list[0]=string2list (buf);
      vfree (buf);
    }
  else if ( file_exists (string[0]))
    {
      list=read_group (string[0]);
    }
  else
    {
      fprintf (stderr, "\nERROR: file <%s> does not exist [FATAL:%s]\n",string[0], PROGRAM);
      exit (EXIT_FAILURE);
    }  

  
  a=0;
  while (list[a])
    {
      int i, b;
      n=atoi (list[a][0]);
      for (b=0; b<A->nseq; b++)nlist[b]=0;
      for (b=2; b<n; b++)
      {
        i=name_is_in_list (list[a][b], A->name, A->nseq, MAXNAMES);
        nlist[i]=1;
      }
      T=collapse_sub_tree ( T,A->nseq,nlist,list[a][1]);
      free_char (list[a], -1);
      a++;
    }
  vfree (list);
  return T;
}

NT_node collapse_sub_tree ( NT_node T,int nseq, int *list, char *new_name)
{
  if (!T) return T;
  else
    {
      int a=0;
 

      while (a<nseq && list[a]==T->lseq2[a]){a++;} 
      if (a==nseq)
      {
        sprintf ( T->name, "%s", new_name);
        T->leaf=T->isseq=1;
        T->left=T->right=NULL;
        return T;
      }
      else
      {
       collapse_sub_tree (T->right, nseq, list, new_name);
       collapse_sub_tree (T->left, nseq, list, new_name);
       return T;
      }
    }
}
  
/*********************************************************************/
/*                                                                   */
/*                                   tree pruning                    */
/*                                                                   */
/*                                                                   */
/*********************************************************************/
NT_node main_prune_tree ( NT_node T, Sequence *S)
{
  prune_tree ( T, S);
  return T;
}

int prune_tree ( NT_node T, Sequence *S)
{
  int nL, nR;
  NT_node C, P;
  
 

  if ( !T)return 0;
  else if ( T->leaf ==1)
    {
      P=T->parent;
      if (name_is_in_list (T->name,S->name, S->nseq, 100)!=-1)return 1;
      else
      {     
      if (P->left==T)  P->left=NULL;
      if (P->right==T) P->right=NULL;
      return 0;
      }
    }
  else
    {
     
      nL=prune_tree (T->left, S);
      nR=prune_tree (T->right, S);
      P=T->parent;
      
      if ( nL> 0 && nR> 0);
      else if (nL==0 && nR==0)/*prune out the node*/
      {
        if ( P==NULL){T->left=NULL;T->right=NULL;}
        else if (P->left ==T)P->left=NULL;
        else if (P->right==T)P->right=NULL;
      }
      else if ( nR==0 || nL==0)/*shunt the current node*/
      {
        
        C=(nR==0)?T->left:T->right;
        
        if      ( P==NULL && nR==0)T->right=NULL;
        else if ( P==NULL && nL==0)T->left=NULL;
        else if ( T==P->left) {P->left=C;}
        else if ( T==P->right){P->right=C;}
        
        if (P)
          {
            C->parent=P;
            C->dist+=T->dist;
          }
      }
    }
  return nL+nR;
}
/*********************************************************************/
/*                                                                   */
/*                                   tree comparison                 */
/*                                                                   */
/*                                                                   */
/*********************************************************************/
NT_node main_compare_trees ( NT_node T1, NT_node T2, FILE *fp)
{
  Sequence *S1, *S2, *S;
  int n;
  float t1, t2, t, w1, w2,w3, w4, w, l1, l;
  
  S1=tree2seq(T1, NULL);
  S2=tree2seq(T2, NULL);  
  S=trim_seq ( S1, S2);

  prune_tree (T1, S);  
  T1=unroot_tree (T1);
  T1=recode_tree(T1, S);
  
  
  prune_tree (T2, S);  
  T2=unroot_tree (T2);
  T2=recode_tree(T2, S);
  
  
  if (getenv4debug("DEBUG_TREE_COMPARE"))fprintf ( stderr, "\ndisplay 1");
  if (getenv4debug("DEBUG_TREE_COMPARE"))display_tree ( T1, S->nseq, stderr);
  if (getenv4debug("DEBUG_TREE_COMPARE"))fprintf ( stderr, "\ndisplay 2");
  if (getenv4debug("DEBUG_TREE_COMPARE"))display_tree ( T2, S->nseq, stderr);
  
  n=tree2nnode(T1);
  t1=compare_trees ( T1, T2, S->nseq, TOPOLOGY);
  t2=compare_trees ( T1, T1, S->nseq, TOPOLOGY);
  t=(t1*100)/t2;

  w1=compare_trees ( T1, T2, S->nseq, WEIGHTED);
  w2=compare_trees ( T1, T1, S->nseq, WEIGHTED);
  w3=compare_trees ( T2, T1, S->nseq, WEIGHTED);
  w4=compare_trees ( T2, T2, S->nseq, WEIGHTED);
  w=((w1+w3)*100)/(w2+w4);
  
  l1=compare_trees ( T1, T2, S->nseq, LENGTH);
  l=100-l1/t2;

  

  
  fprintf ( fp, "\n#tree_cpm|T: %.f W: %.2f L: %.2f", t, w, l);
  fprintf ( fp, "\n#tree_cpm|%d Nodes in T1 with %d Sequences",n, S->nseq);
  fprintf ( fp, "\n#tree_cmp|T: ratio of identical nodes");
  fprintf ( fp, "\n#tree_cmp|W: ratio of identical nodes weighted with the min Nseq below node");
  fprintf ( fp, "\n#tree_cmp|L: average branch length similarity\n");
  
  compare_trees ( T1, T2, S->nseq, RECODE);
  
  free_sequence (S, -1);
  free_sequence (S1, -1);
  free_sequence (S2, -1);
  
  return T1;
}  
  
float compare_trees ( NT_node T1, NT_node T2, int nseq,int  mode)
{
  /*search each branch of T1 in T2*/
  float n=0;
  
  
  if ( !T1 || !T2)return 0;
  
  if (getenv4debug("DEBUG_TREE_COMPARE"))display_node (T1, "\nNODE ", nseq);
  n+=search_node ( T1, T2, nseq, mode);
  
  n+=compare_trees ( T1->left, T2, nseq, mode);
  n+=compare_trees ( T1->right, T2, nseq, mode);
  n+=compare_trees ( T1->bot, T2, nseq, mode);
  return n;
}

float search_node ( NT_node B, NT_node T, int nseq, int mode)
{
  int n=0;
  if ( !B || !T) return -1;
  if (getenv4debug("DEBUG_TREE_COMPARE"))display_node ( T, "\n\t", nseq);
  
  n=compare_node ( B->lseq2, T->lseq2, nseq );
  
  if ( n==1) 
    {
      if (getenv4debug("DEBUG_TREE_COMPARE"))fprintf ( stderr, "[1][%d]", (int)evaluate_node_similarity ( B, T, nseq, mode));
      if (mode==RECODE)B->dist=B->leaf;
      return evaluate_node_similarity ( B, T, nseq, mode);
    }
  else if ( n==-1)
    {
      if (getenv4debug("DEBUG_TREE_COMPARE"))fprintf ( stderr, "[-1]"); 
      if (mode==RECODE)B->dist=-B->leaf;
      return 0;
    }
  else
    {
      if (getenv4debug("DEBUG_TREE_COMPARE"))fprintf ( stderr, "[0]");
      n=search_node ( B, T->left, nseq, mode);   
      if ( n>0) return n;
      n=search_node ( B, T->right, nseq, mode);
      if ( n>0) return n;
      n=search_node ( B, T->bot, nseq, mode);
      if ( n>0) return n;
    }
  return n;
}

float evaluate_node_similarity ( NT_node B, NT_node T, int nseq, int mode)
{
int a, c;
 
 if ( mode==TOPOLOGY || mode ==RECODE)
   {
     for ( a=0; a< nseq; a++)
       if ( B->lseq2[a]!=T->lseq2[a]) return 0;
     return 1;
   }
 else if ( mode == WEIGHTED)
   {
     for (c=0, a=0; a< nseq; a++)
       {
       if ( B->lseq2[a]!=T->lseq2[a]) return 0; 
       else c+=B->lseq2[a];
       }      
     return (float)(MIN(c,nseq));
   }
 else if ( mode == LENGTH )
   {
     float d1, d2;
     
     for (c=0, a=0; a< nseq; a++)
       {
       if ( B->lseq2[a]!=T->lseq2[a]) return 0; 
       }
     d1=FABS((B->dist-T->dist));
     d2=MAX(B->dist, T->dist);
     return (d2>0)?(d1*100)/d2:0;
   }
 else
   {
     return 0;
   }
}
int compare_node ( int *b1, int *b2, int nseq)
{
  int n1, n2;
  
  n1=compare_node1 ( b1, b2, nseq);
  /*fprintf ( stderr, "[%d]", n1);*/
  if ( n1==1) return 1;
  
  n2=compare_node2 ( b1, b2, nseq);
  /* fprintf ( stderr, "[%d]", n2);*/
  if ( n2==1)return 1;
  else if ( n2==-1 && n1==-1) return -1;
  else return 0;
}
int compare_node1 ( int *b1, int *b2, int n)
{
  int a;
  int l1, l2;
  int r=1;
  for ( a=0; a< n; a++)
    {
      l1=b1[a];
      l2=b2[a];
      if ( l1==1 && l2==0) return -1;
      if ( l1!=l2)r=0;
    }
  return r;
}
int compare_node2 ( int *b1, int *b2, int n)
{
  int a;
  int l1, l2;
  int r=1;
  
  for ( a=0; a< n; a++)
    {
      l1=1-b1[a];
      l2=b2[a];
      if ( l1==1 && l2==0) return -1;
      if ( l1!=l2) r=0;
    }
  return r;
}
void display_node (NT_node N, char *string,int nseq)
{
  int a;
  fprintf ( stderr, "%s", string);
  for (a=0; a< nseq; a++)fprintf ( stderr, "%d", N->lseq2[a]);
  
}

NT_node recode_tree ( NT_node T, Sequence *S)
{
  int a;
  /*
    ->lseq2[N] is a coded version of the content of each node
    ->lseq2[N] contains 0 if sequence n is in the node 
    ->lseq2[N] contains 1 if the sequendce is under the considered node
    This function recodes lseq2 so that the index n is given by the sequence order in the datastructure S
  */
  
  if ( !T)return NULL;
  
  
  vfree ( T->lseq2);T->lseq2=vcalloc (S->nseq+1, sizeof (int));
  
  if ( T->leaf==1)
    {
      if((a=name_is_in_list (T->name,S->name, S->nseq, 100))!=-1)
      {
        T->lseq[0]=a;
        T->lseq2[a]=1;
      }
      else
      {
        fprintf (stderr, "\nPB with sequence %s", T->name);
      }
    }
  else
    {
      T->left=recode_tree (T->left, S);
      T->right=recode_tree (T->right,S);
      T->bot=recode_tree (T->bot, S);
      T->nseq=0;
      if ( T->left)
      {
        for ( a=0; a<(T->left)->nseq; a++, T->nseq++)
          T->lseq[T->nseq]=(T->left)->lseq[a];
        for ( a=0; a<S->nseq; a++)
          if ( (T->left)->lseq2[a]) T->lseq2[a]=1;
        
      }
      if ( T->right)
      {
        for ( a=0; a<(T->right)->nseq; a++, T->nseq++)
          T->lseq[T->nseq]=(T->right)->lseq[a];
        for ( a=0; a<S->nseq; a++)
          if ( (T->right)->lseq2[a]) T->lseq2[a]=1;
      }
      if ( T->bot)
      {
        for ( a=0; a<(T->bot)->nseq; a++, T->nseq++)
          T->lseq[T->nseq]=(T->bot)->lseq[a];
        for ( a=0; a<S->nseq; a++)
          if ( (T->bot)->lseq2[a]) T->lseq2[a]=1;
      }
    }
  return T;
}


/*********************************************************************/
/*                                                                   */
/*                                   FJ_tree Computation             */
/*                                                                   */
/*                                                                   */
/*********************************************************************/
static NT_node compute_fj_tree (NT_node T, Alignment *A, int limit, char *mode);
static NT_node compute_cw_tree (Alignment *A);
static NT_node compute_std_tree (Alignment *A, int n, char **arg);
NT_node tree_compute ( Alignment *A, int n, char ** arg_list)
{
  if (n==0 || strm (arg_list[0], "cw"))
    {
      return compute_cw_tree (A);
    }
  else if ( strm (arg_list[0], "fj"))
    {
      return compute_fj_tree ( NULL, A, atoi(arg_list[0]), (n>=2)?arg_list[1]:NULL);
    }
  
  else if ( ( strm (arg_list[0], "nj")))
    {
      return compute_std_tree (A, n, arg_list);
    }
  else
    return compute_std_tree (A, n, arg_list);
}
NT_node compute_std_tree (Alignment *A, int n, char **arg_list)
{
  int **s, i;  
  NT_node T;
  char *tree_name;
  
  tree_name =vtmpnam (NULL);
    
  if ( name_is_in_list ("aln",arg_list,n, 3)!=-1)
    {
      s=get_sim_aln_array ( A, "idmat");
    }
  else if ((i=name_is_in_list ("ktup", arg_list,n,4))!=-1)
    {
      ungap_array (A->seq_al, A->nseq);
      s=get_sim_aln_array ( A,arg_list[i] );
    }
  
  s=sim_array2dist_array(s, 100);
  int_dist2nj_tree (s, A->name, A->nseq, tree_name);
  T=main_read_tree(tree_name);

  if ( name_is_in_list ("dpa", arg_list,n, 3)!=-1)
    {
      s=dist_array2sim_array(s, 100);
      T=code_dpa_tree (T,s);
    }
  return T;
}
        
  
  
NT_node compute_cw_tree (Alignment *A)
{
  char *tmp1, *tmp2, tmp3[1000];
  char command[1000];
  
  tmp1=vtmpnam (NULL);
  tmp2=vtmpnam (NULL);
  sprintf ( tmp3, "%s.ph", tmp1);
  output_clustal_aln (tmp1, A);
  sprintf ( command, "clustalw -infile=%s -tree %s", tmp1, TO_NULL_DEVICE);
  my_system ( command);
  sprintf ( command, "mv %s %s", tmp3, tmp2);
  my_system ( command);
  return main_read_tree(tmp2);
}

NT_node compute_fj_tree (NT_node T, Alignment *A, int limit, char *mode)
{
  static int in_fj_tree;
  if (!in_fj_tree)fprintf ( stderr, "\nComputation of an NJ tree using conserved positions\n");
  
  in_fj_tree++;
  if (T && T->leaf<=2);
  else
    {
      T=aln2fj_tree(T,A,limit, mode);
      T->right=compute_fj_tree ( T->right, A, limit, mode);
      T->left=compute_fj_tree  ( T->left, A, limit, mode);
    }
  in_fj_tree--;
  return T;
}

NT_node aln2fj_tree(NT_node T, Alignment *A, int limit, char *mode)
{
  NT_node NT;
  Sequence *S=NULL;
  Alignment *subA;
  
  if (T)S=tree2seq (T,NULL);
 

  subA=extract_sub_aln2 (A,(T)?S->nseq:A->nseq, (T)?S->name:A->name);
  subA=filter_aln4tree (subA,limit,mode);
  
  
  while ( subA->len_aln<1)
    {
      free_aln (subA);
      subA=extract_sub_aln2 (A,(T)?S->nseq:A->nseq, (T)?S->name:A->name);
      subA=filter_aln4tree (subA,--limit, mode);
    }

  NT=aln2tree (subA);

  NT=realloc_tree (NT,A->nseq); 
  fprintf ( stderr, "Limit:%d Columns: %4d Left: %4d Right %4d BL:%4.2f\n",limit, subA->len_aln, (NT->right)->leaf,(NT->left)->leaf, (NT->left)->dist+(NT->right)->dist);
  
  if ( T)
    {
      NT->dist=T->dist;
      NT->parent=T->parent;
    }
  free_tree(T);
  free_aln (subA);
  free_sequence (S, -1);
  return NT;
}  
  
Alignment * filter_aln4tree (Alignment *A, int n, char *mode)
{
  char *aln_file;
  char *score_file;
  char *saln_file;
  char command[1000];
  
  aln_file=vtmpnam(NULL);
  score_file=vtmpnam(NULL);
  saln_file=vtmpnam(NULL);
  
  output_clustal_aln (aln_file, A);
 
  sprintf ( command, "%s -in %s -action +evaluate %s -output clustalw > %s", getenv ("SEQ_REFORMAT"),aln_file,(mode)?mode:"categories", score_file);
  my_system ( command);
  
  sprintf ( command, "%s -in %s -struc_in %s -struc_in_f number_aln -action +use_cons +keep '[%d-8]' +rm_gap -output clustalw > %s", getenv("SEQ_REFORMAT"), aln_file, score_file,n, saln_file);
  my_system ( command);
  
  system ( command);
  free_aln (A);
  
  A=main_read_aln ( saln_file, NULL);
  return A;
}


/*********************************************************************/
/*                                                                   */
/*                                   Tree Filters and MAnipulation   */
/*                                                                   */
/*                                                                   */
/*********************************************************************/
NT_node aln2tree (Alignment *A)
{
  NT_node **T=NULL;
  

  T=make_nj_tree (A, NULL, 0, 0, A->seq_al, A->name, A->nseq, NULL, NULL);
  tree2nleaf (T[3][0]);
  
  return T[3][0];
}
NT_node realloc_tree ( NT_node R, int n)
{
  if ( !R)return R;
  if ( !R->leaf) 
    {
      R->right=realloc_tree (R->right,n);
      R->left=realloc_tree (R->left,n);
      R->bot=realloc_tree (R->bot,n);
    }
  R->lseq=vrealloc (R->lseq, n*sizeof (int));
  R->lseq2=vrealloc (R->lseq2, n*sizeof (int));
  return R;
}

NT_node reset_boot_tree ( NT_node R, int n)
{
  if ( !R)return R;
  if ( !R->leaf) 
    {
      
      R->right=reset_boot_tree (R->right,n);
      R->left=reset_boot_tree (R->left,n);
      R->bot=reset_boot_tree (R->bot,n);
    }
  R->bootstrap=(float)n;
          
  return R;
}

NT_node reset_dist_tree ( NT_node R, float n)
{
  if ( !R)return R;
  if ( !R->leaf) 
    {
      
      R->right=reset_dist_tree (R->right,n);
      R->left=reset_dist_tree (R->left,n);
      R->bot=reset_dist_tree (R->bot,n);
    }
  if (R->parent && !(R->parent)->parent && !(R->parent)->bot)R->dist=n/2;
  else R->dist=n;

  return R;
}

NT_node free_tree ( NT_node R)
{
  if ( !R)return R;
  if ( R->leaf!=1) 
    {
      R->right=free_tree (R->right);
      R->left=free_tree (R->left);
      R->bot=free_tree (R->bot);
    }
  free_tree_node (R);
  return R;
}

NT_node free_tree_node ( NT_node R)
{
  vfree ( R->name);
  vfree ( R->lseq); vfree ( R->lseq2);
  vfree (R);
  return NULL;
}

NT_node   rename_seq_in_tree ( NT_node R, char ***list)
{
  if ( !R || !list) return R;

  if ( R->leaf!=1) 
    {
      R->right=rename_seq_in_tree (R->right, list);
      R->left=rename_seq_in_tree (R->left, list);
      R->bot=rename_seq_in_tree (R->bot, list);
    }
  else
    {
      int n=0;
      while ( list[n][0][0])
      {
        if ( strm (list[n][0], R->name))sprintf (R->name, "%s",list[n][1]);
        n++;
      }
    }
  return R;
}
Sequence * tree2seq    (NT_node R, Sequence *S)
{

  if ( !R)return S;
  if ( !S)
    {
      S=declare_sequence (10, 10, tree2nseq (R));
      S->nseq=0;
    }
  
  if (R->leaf==1)
    {
      sprintf ( S->name[S->nseq++], "%s", R->name);
    }
  else
    {
      S=tree2seq (R->left, S);
      S=tree2seq (R->right, S);
      S=tree2seq (R->bot, S);
    }
  return S;
}

NT_node balance_tree (NT_node T)
{
  static int **list;
  NT_node NL[3];
  
  if ( !T) return T;
  else if ( T->leaf<=2)return T;
  else    
    {
      if (!list)list=declare_int (3, 2);
      
      NL[0]=T->left;
      NL[1]=T->right;
      NL[2]=T->bot;
      
      list[0][0]=(T->left)?(T->left)->leaf:0;
      list[0][1]=0;
      list[1][0]=(T->right)?(T->right)->leaf:0;
      list[1][1]=1;
      list[2][0]=(T->bot)?(T->bot)->leaf:0;
      list[2][1]=2;
      
      sort_int (list,2,0,0,2);
      
      T->left=NL[list[2][1]];
      T->right=NL[list[1][1]];
      T->bot=NL[list[0][1]];
      
      T->left=balance_tree (T->left);
      T->right=balance_tree (T->right);
      T->bot=balance_tree (T->bot);
      return T;   
    }
}
FILE * display_tree (NT_node R, int nseq, FILE *fp)
{
  int a;
  
  if ( !R);
  else 
    {
      /*
      if ( R->nseq==1)fprintf (stderr,"\n[%s] ", R->name);
      else fprintf ( stderr, "\n[%d Node] ",R->nseq);
      for ( a=0; a< R->nseq; a++) fprintf ( stderr, "[%d]", R->lseq[a]);
      */
      fprintf (fp, "\n");
      for ( a=0; a< nseq; a++)fprintf (fp, "%d", R->lseq2[a]);
      if (R->leaf==1) fprintf (fp, " %s", R->name);
      fprintf (fp, " :%.4f", R->dist);
      fp=display_tree (R->left, nseq, fp);
      fp=display_tree (R->right, nseq, fp);
      fp=display_tree (R->bot, nseq, fp);
    }
  return fp;
}

int tree2nnode ( NT_node R)
{
  int n;
  if ( !R)n=0;
  else if ( R->leaf==1){R->node=1;n=1;}
  else 
    {
      n=1;
      n+=tree2nnode (R->right);
      n+=tree2nnode (R->left);
      n+=tree2nnode (R->bot);
      R->node=n;
    }
  return n;
}
int tree2nleaf (NT_node R)
{
  if ( !R)return 0;
  else if (R->leaf==1){return 1;}
  else if (R->right==NULL && R->left==NULL && R->bot==NULL){R->leaf=1; return 1;}
  else
    {
      int n=0;
      n+=tree2nleaf (R->right);
      n+=tree2nleaf (R->left);
      n+=tree2nleaf (R->bot);
      
      R->leaf=n;
      return n;
    }
}


int tree2nseq ( NT_node R)
{
 return tree2nleaf(R);
}

int tree_file2nseq (char *fname)
{
  FILE *fp;
  char *string;
  int p, a, b, c, n;
  
  string=vcalloc (count_n_char_in_file(fname)+1, sizeof (char));
  
  fp=vfopen (fname, "r");
  n=0;
  while ( (c=fgetc(fp))!=EOF){if (c=='(' || c==')' || c==',' || c==';') string[n++]=c;}
  vfclose (fp);string[n]='\0';
  
  for (n=0, p=1; p<strlen (string); p++)
    {
      a=string[p-1];
      b=string[p];

      if ( a=='(' && b==',')n++;
      else if ( a==',' && b==')')n++;
      else if ( a==',' && b==',')n++;
    }
  if ( getenv4debug("DEBUG_TREE"))fprintf (stderr, "\n[DEBUG_TREE:tree_file2nseq]%s", string);
  vfree (string);
  return n;
}
  

void clear_tree ( NT_node R)
{
  if (!R) return;
  
  R->visited=0;
  
  if ( R->leaf==1);
  else
    {
      clear_tree ( R->right);
      clear_tree ( R->left);
      clear_tree ( R->bot);
    }
}
void display_leaf ( NT_node T)
{
  if ( !T)return;
  else if ( T->visited)return;
  else T->visited=1;
  
  if ( T->leaf==1)fprintf ( stderr, " %s", T->name);
  else
    {
      display_leaf ( T->right);
      display_leaf ( T->left);
      display_leaf ( T->bot);
    }
}

NT_node find_longest_branch ( NT_node T, NT_node L)
  {
    
    if ( !L || T->dist>L->dist)
      {
      L=T;
      }

    if ( T->leaf==1)return L;
    else 
      {
      L=find_longest_branch ( T->right, L);
      L=find_longest_branch ( T->left,  L);
      return L;
      }
  }

NT_node reroot_tree ( NT_node TREE, NT_node R)
{
  /*ReRoots the tree between Node T and its parent*/
  NT_node NewRoot;
  NT_node L;
  float dist;
  
  if ( !R->parent)return R;
  TREE=unroot_tree (TREE);
  NewRoot=declare_tree_node(TREE->maxnseq);
  
  L=R->parent;
  
  NewRoot->right=R;
  NewRoot->left=L;
  dist=R->dist/2;

  R->parent=NULL;

  if ( L->left==R)L->left=NULL;
  else if ( L->right==R)L->right=NULL;
  else if ( L->bot==R)L->bot=NULL;
  
  NewRoot->left=straighten_tree ( NewRoot, NewRoot->left, dist);
  NewRoot->right=straighten_tree ( NewRoot, NewRoot->right, dist);

  tree2nleaf(NewRoot);
  return NewRoot;
}
  
  
  
NT_node straighten_tree ( NT_node P, NT_node C, float new_dist)
{
  float dist;
  int n;
  NT_node list[2];
  
  if ( C==NULL)return NULL;


  fprintf ( stderr, "\nN Leaf: %d", C->leaf);
  dist=C->dist;
  C->dist=new_dist;
  n=0;
  list[0]=list[1]=NULL;
  if ( C->right  && C->right!=P ) {fprintf ( stderr, "\nR");list[n++]=C->right;}
  if ( C->left   && C->left!=P  ) {fprintf ( stderr, "\nL");list[n++]=C->left;}
  if ( C->parent && C->parent!=P) {fprintf ( stderr, "\nP");list[n++]=C->parent;}
  if ( C->bot    && C->bot!=P   ) {fprintf ( stderr, "\nB");list[n++]=C->bot;}
  C->right=C->left=C->parent=C->bot=NULL;

  
  C->parent=P;
  C->dist=new_dist;
  
  C->right=list[0];
  C->left= list[1];

  
  if ( n>2) fprintf ( stderr, "\nERROR: too Many Nodes");
 
  
  C->right=straighten_tree (C, C->right, dist);
  C->left=straighten_tree  (C, C->left,  dist);
  return C;
}


NT_node unroot_tree ( NT_node T)
{
  
  if ( !T) return NULL;
  else if ( T->left==NULL && T->left==NULL)return NULL;
  else if ( T->left && T->right && T->bot && T->parent==NULL)return T;
  else if ( T->left==NULL || T->right==NULL)
    {
      NT_node C;
      C=(T->left==NULL)?T->right:T->left;
      C->parent=NULL;
      C->dist=0;
      
      return unroot_tree (C);
    }
  else
    {
      float dist;
      clear_tree (T);
      dist=(T->right)->dist+(T->left)->dist;
      (T->right)->dist=dist;
      (T->right)->parent=T->left;
      (T->left)->bot=(T->right);
      (T->left)->dist=0;
      (T->left)->parent=NULL;
    }
  return T->left;
}



FILE * print_tree ( NT_node T, char *format,FILE *fp)
{
  Sequence *S;
  
  tree2nleaf(T);
  S=tree2seq(T, NULL);
  recode_tree (T, S);
  
  if ( strm (format, "binary"))
       fp=display_tree ( T,S->nseq, fp);
  else if ( strm (format, "newick"))
    {
      T=balance_tree (T);
      fp=rec_print_tree (T, fp);
      fprintf ( fp, ";\n");
    }
  else
    {
      fprintf ( stderr, "\nERROR: %s is an unknown tree format [FATAL:%s]\n", format, PROGRAM);
      exit (EXIT_FAILURE);
    }
  return fp;
}
FILE * rec_print_tree ( NT_node T, FILE *fp)
{
  
  if (T==NULL || T->visited==1)
    return fp;
  else 
    T->visited=1;
  
  if ( T->isseq)
    {
      fprintf ( fp, " %s:%.5f",T->name, T->dist);
    }
  else
    {
      fprintf ( fp, "(");fp=rec_print_tree ( T->left, fp);
      fprintf ( fp, ",");fp=rec_print_tree ( T->right, fp);
      if ( T->bot){fprintf ( fp, ",");fp=rec_print_tree ( T->bot, fp);}
      fprintf ( fp, ")");
      if (T->parent || T->dist)
      {
        if ( T->bootstrap!=0)fprintf (fp, " %.2f", T->bootstrap);
        fprintf (fp, ":%.5f", T->dist);
      }
    }
  return fp;
}

/*********************************************************************/
/*                                                                   */
/*                                  Tree Functions                   */
/*                                                                   */
/*                                                                   */
/*********************************************************************/

int ** make_sub_tree_list ( NT_node **T, int nseq, int n_node)
        {


/*This function produces a list of all the sub trees*/
       

/*      /A */
/*    -* */
/*      \      /B */
/*       \    / */
/*        ---* */
/*                  \ */
/*                   *--C */
/*                    \ */
/*                     \D */
      
/*    Contains 4 i_nodes */
/*           8   nodes (internal nodes +leaves) */
/*           8 sub trees: */
/*           ABCD */
/*                 1111 */
/*           0111 */
/*           1000 */
/*           0100 */
/*           0011 */
/*           0001 */
/*           0010 */
       
      int **sub_tree_list;
      int a, n=0;
      
      
      if (T)
           {
             sub_tree_list=declare_int ( (n_node), nseq);
             make_all_sub_tree_list (T[3][0],sub_tree_list, &n);

           }
      else
           {
             sub_tree_list=declare_int (nseq, nseq);
             for ( a=0; a< nseq; a++)sub_tree_list[a][a]=1;
           }
      
      return sub_tree_list;
      }

void make_all_sub_tree_list ( NT_node N, int **list, int *n)
        {
      make_one_sub_tree_list (N, list[n[0]++]);
      if (N->leaf!=1)
          {
          make_all_sub_tree_list (N->left , list, n);
          make_all_sub_tree_list (N->right, list, n);
          }
      return;
      }

void make_one_sub_tree_list ( NT_node T,int *list)
        {
      if (T->leaf==1)
          {
          
          list[T->seq]=1;
          }
      else
          {
          make_one_sub_tree_list(T->left , list);
          make_one_sub_tree_list(T->right, list);
          }
      return;
      }


NT_node main_read_tree(char *treefile)
{
  int tot_node=0;
  NT_node **T;
  
  T=read_tree ( treefile, &tot_node,tree_file2nseq (treefile),NULL);
  
  return T[3][0];
}


NT_node** read_tree(char *treefile, int *tot_node,int nseq, char **seq_names)
      {

          /*The Tree Root is in the TREE[3][0]...*/
          /*TREE[0][ntot]--> pointer to each node and leave*/
      char ch;
      int  a,b;

      FILE *fp;
      int nseq_read = 0;
      int nnodes = 0;/*Number of Internal Nodes*/
      int ntotal = 0;/*Number of Internal Nodes + Number of Leaves*/
      int flag;
      int c_seq;
      NT_node **lu_ptr;
      NT_node seq_tree, root,p;
      

      tot_nseq=nseq;
      rooted_tree=distance_tree=TRUE;
      
      fp = vfopen(treefile, "r");
      fp=skip_space(fp);
      ch = (char)getc(fp);
      if (ch != '(')
            {
            fprintf(stderr, "Error: Wrong format in tree file %s\n", treefile);
                  exit (1);
            }
      rewind(fp);
  
      
      lu_ptr=(NT_node **)vcalloc(4,sizeof(NT_node*));
      lu_ptr[0] = (NT_node *)vcalloc(10*nseq,sizeof(NT_node));
      lu_ptr[1] = (NT_node *)vcalloc(10*nseq,sizeof(NT_node));
      lu_ptr[2] = (NT_node *)vcalloc(10*nseq,sizeof(NT_node));
      lu_ptr[3] =(NT_node *) vcalloc(1,sizeof(NT_node));
      
      seq_tree =(NT_node) declare_tree_node(nseq);
  
      set_info(seq_tree, NULL, 0, "  ", 0.0, 0);
      
      fp=create_tree(seq_tree,NULL,&nseq_read, &ntotal, &nnodes, lu_ptr, fp);
      fclose (fp);
            
      
      if (nseq != tot_nseq)
            {
            fprintf(stderr," Error: tree not compatible with alignment (%d sequences in alignment and %d in tree\n", nseq,nseq_read);
            exit (1);
            }
      if (distance_tree == FALSE)
            {
            if (rooted_tree == FALSE)
                  {
                        fprintf(stderr,"Error: input tree is unrooted and has no distances, cannot align sequences\n");
                        exit (1);
                  }
            }
        if (rooted_tree == FALSE)
            {

              root = reroot(seq_tree, nseq,ntotal,nnodes, lu_ptr);
              lu_ptr[1][nnodes++]=lu_ptr[0][ntotal++]=root;
              
            }
      else
            {
            root = seq_tree;
            }     
      lu_ptr[3][0]=root;
      tot_node[0]=nnodes;
      
      
      
      for ( a=0; a< ntotal; a++)
            {
            (lu_ptr[0][a])->isseq=(lu_ptr[0][a])->leaf;
            (lu_ptr[0][a])->dp=(lu_ptr[0][a])->dist;
            }
            
      
      for ( a=0; a< nseq; a++)
            {
            if (!seq_names)
              {
                flag=1;
                (lu_ptr[2][a])->order=(lu_ptr[2][a])->seq=a;
              }
            else
              {
                for ( flag=0,b=0; b<nseq; b++)
                  {                  
                  if ( strncmp ( (lu_ptr[2][a])->name, seq_names[b], MAXNAMES)==0)
                    {
                      flag=1;
                      
                      (lu_ptr[2][a])->order=(lu_ptr[2][a])->seq=b;
                        /*vfree ( (lu_ptr[2][a])->name);*/
                      sprintf ((lu_ptr[2][a])->name, "%s", seq_names[b]);
                    }
                  }
              }
            if ( flag==0  && (lu_ptr[0][a])->leaf==1)
                  {
                  fprintf ( stderr, "\n%s* not in tree",(lu_ptr[2][a])->name);
                  for ( a=0; a< ntotal; a++)
                    {
                      fprintf ( stderr, "\n%d %s",(lu_ptr[2][a])->leaf, (lu_ptr[2][a])->name);
                    } 
                  } 
            }
      
      for ( a=0; a< nseq; a++)
            {
            p=lu_ptr[2][a];
            c_seq=p->seq;
            
            while ( p!=NULL)
                  {
                  p->lseq[p->nseq]=c_seq;
                  p->nseq++;
                  p=p->parent;
                  }
            }
      
      
      return lu_ptr;    
      }


FILE * create_tree(NT_node ptree, NT_node parent,int *nseq,int  *ntotal,int  *nnodes,NT_node **lu, FILE *fp)
      {

      int i, type;
      float dist=0;
      float bootstrap=0;
      char *name;
      int ch;
      

      name=vcalloc ( MAXNAMES+1, sizeof (char));
      sprintf ( name, "   ");
      fp=skip_space(fp);
      ch = (char)getc(fp);
      
      if (ch == '(')
            {
            type = NODE;
                  name[0] = '\0';
                  lu[0][ntotal[0]] = lu[1][nnodes[0]] = ptree;
                  ntotal[0]++;
                  nnodes[0]++;
                  create_tree_node(ptree, parent);
                  fp=create_tree(ptree->left, ptree, nseq,ntotal,nnodes,lu,fp);
            ch = (char)getc(fp);
                if ( ch == ',')
                        {
                  fp=create_tree(ptree->right, ptree,nseq,ntotal,nnodes,lu,fp);
                  ch = (char)getc(fp);
                  if ( ch == ',')
                              {
                              
                              ptree = insert_tree_node(ptree);
                              lu[0][ntotal[0]] = lu[1][nnodes[0]] = ptree;
                              ntotal[0]++;
                              nnodes[0]++;
                              fp=create_tree(ptree->right, ptree,nseq,ntotal,nnodes,lu,fp);
                              rooted_tree = FALSE;
                        if ( getenv4debug ( "DEBUG_TREE")){fprintf ( stderr, "\n[DEBUG_TREE:create_tree] Unrooted Tree");}
                              }
                        }

                  fp=skip_space(fp);
                  ch = (char)getc(fp);
            }
      else
            {
            type=LEAF;
                  lu[0][ntotal[0]] = lu[2][nseq[0]] = ptree;
                  ntotal[0]++;
                  nseq[0]++;
                  name[0] = ch;
                  i=1;
                  ch = (char)getc(fp);
                  while ((ch != ':') && (ch != ',') && (ch != ')'))
                  {
                  if (i < MAXNAMES) name[i++] = ch;
                  ch = (char)getc(fp);
                  }
                  name[i] = '\0';
            
                  if ( i>=(MAXNAMES+1)){fprintf (stderr, "\nName is too long");exit (1);}
                  if (ch != ':' && !isdigit(ch))
                  {
                    /*distance_tree = FALSE*/;
                  }
            }
      if (ch == ':')
            {
                  fp=skip_space(fp);
                  fscanf(fp,"%f",&dist);
                  fp=skip_space(fp);
            bootstrap=0;
                  }
      /*Tree with Bootstrap information*/
      else if (isdigit (ch))
        {
          ungetc(ch,fp);
          fscanf(fp,"%f",&bootstrap);
          if ( fscanf(fp,":%f",&dist)==1);
          else dist=0;
          fp=skip_space(fp);
        }
      else
        {
          ungetc ( ch, fp);
          skip_space(fp);
        }
        set_info(ptree, parent, type, name, dist, bootstrap);
      
      
      vfree (name);
      return fp;
      }

NT_node declare_tree_node (int nseq)
      {
      NT_node p;
      
      p= (NT_node)vcalloc (1, sizeof ( Treenode)); 
      p->left = NULL;
      p->right = NULL;
      p->parent = NULL;
      p->dist = 0.0;
      p->leaf = 0;
      p->order = 0;
      p->maxnseq=nseq;
      p->name=(char*)vcalloc (MAXNAMES+1,sizeof (char));
      p->name[0]='\0';
      p->lseq=(int*)vcalloc ( nseq, sizeof (int));
      return p;
      
      }

void set_info(NT_node p, NT_node parent, int pleaf, char *pname, float pdist, float bootstrap)
      {
      p->parent = parent;
      p->leaf = pleaf;
      p->dist = pdist;
      p->bootstrap=bootstrap;
      p->order = 0;
      
      sprintf (p->name, "%s", pname);
      if (pleaf ==1)
            {
            p->left = NULL;
            p->right = NULL;
            }

   }
NT_node insert_tree_node(NT_node pptr)
      {

      NT_node newnode;

      newnode = declare_tree_node( pptr->maxnseq);
      create_tree_node(newnode, pptr->parent);

      newnode->left = pptr;
      pptr->parent = newnode;

      set_info(newnode, pptr->parent, 0, "", 0.0, 0);

      return(newnode);
      }

void create_tree_node(NT_node pptr, NT_node parent)
      {
      pptr->parent = parent;
      pptr->left =declare_tree_node(pptr->maxnseq) ;
      (pptr->left)->parent=pptr;
      
      pptr->right =declare_tree_node(pptr->maxnseq) ;    
      (pptr->right)->parent=pptr;
      }     

FILE * skip_space(FILE *fp)
      {
      int   c;
  
      do
      c = getc(fp);
      while(isspace(c));
      if ( c==EOF)
            {
            fprintf ( stderr, "\nEOF");
            exit (1);
            }
      ungetc(c, fp);
      return fp;
      }


NT_node reroot(NT_node ptree, int nseq, int ntotal, int nnodes, NT_node **lu)
      {
      NT_node p, rootnode, rootptr;
      float   diff, mindiff, mindepth = 1.0, maxdist;
      int   i;
      int first = TRUE;
      
      
      
      rootptr = ptree;
      
      
      for (i=0; i<ntotal; i++)
            {
            p = lu[0][i];
            if (p->parent == NULL)
                  diff = calc_root_mean(p, &maxdist, nseq, lu);
            else
                  diff = calc_mean(p, &maxdist, nseq, lu);

            if ((diff == 0) || ((diff > 0) && (diff < 2 * p->dist)))
                   {
                        if ((maxdist < mindepth) || (first == TRUE))
                        {
                              first = FALSE;
                              rootptr = p;
                              mindepth = maxdist;
                              mindiff = diff;
                        }
                  }

            }
      
      if (rootptr == ptree)
            {
            mindiff = rootptr->left->dist + rootptr->right->dist;
            rootptr = rootptr->right;
            }
      
      rootnode = insert_root(rootptr, mindiff);
      
      diff = calc_root_mean(rootnode, &maxdist, nseq, lu);

      return(rootnode);
      }


float calc_root_mean(NT_node root, float *maxdist, int nseq, NT_node **lu)
      {
      float dist , lsum = 0.0, rsum = 0.0, lmean,rmean,diff;
      NT_node p;
      int i;
      int nl, nr;
      int direction;
      
   
      dist = (*maxdist) = 0;
      nl = nr = 0;
      for (i=0; i< nseq; i++)
        {
          p = lu[2][i];
          dist = 0.0;
          while (p->parent != root)
            {
                  dist += p->dist;
                  p = p->parent;
            }
         if (p == root->left) direction = LEFT;
         else direction = RIGHT;
         dist += p->dist;

         if (direction == LEFT)
            {
            lsum += dist;
                  nl++;
            }
         else
            {
                  rsum += dist;
                  nr++;
            }
        
        if (dist > (*maxdist)) *maxdist = dist;
      }

      lmean = lsum / nl;
      rmean = rsum / nr;

      diff = lmean - rmean;
      return(diff);
      }

float calc_mean(NT_node nptr, float *maxdist, int nseq,NT_node **lu)
      {
      float dist , lsum = 0.0, rsum = 0.0, lmean,rmean,diff;
      NT_node p, *path2root;
      float *dist2node;
      int depth = 0, i,j , n;
      int nl , nr;
      int direction, found;
      
      
      path2root = (NT_node *)vcalloc(nseq,sizeof(Treenode));
      dist2node = (float *)vcalloc(nseq,sizeof(float));

      depth = (*maxdist) = dist = 0;
      nl = nr = 0;
      p = nptr;
      while (p != NULL)
            {
            path2root[depth] = p;
            dist += p->dist;
            dist2node[depth] = dist;
            p = p->parent;
            depth++;
            }
 
/***************************************************************************
   *nl = *nr = 0;
   for each leaf, determine whether the leaf is left or right of the node.
   (RIGHT = descendant, LEFT = not descendant)
****************************************************************************/
      for (i=0; i< nseq; i++)
            {
                  p = lu[2][i];
                  if (p == nptr)
                  {
                        direction = RIGHT;
                        dist = 0.0;
                  }
                  else
                  {
                  direction = LEFT;
                  dist = 0.0;
                  
                  found = FALSE;
                  n = 0;
                  while ((found == FALSE) && (p->parent != NULL))
                        {
                              for (j=0; j< depth; j++)
                              if (p->parent == path2root[j])
                                          { 
                                          found = TRUE;
                                          n = j;
                                          }
                              dist += p->dist;
                              p = p->parent;
                        }
         
                  if (p == nptr) direction = RIGHT;
        
                  }
            if (direction == LEFT)
                  {
                        lsum += dist;
                        lsum += dist2node[n-1];
                        nl++;
                  }
            else
                  {
                        rsum += dist;
                        nr++;
                  }

            if (dist > (*maxdist)) *maxdist = dist;
            }

   vfree(dist2node);
   vfree(path2root);

   

   if ( nl==0 || nr==0)exit (1);    
   lmean = lsum / nl;
   rmean = rsum / nr;
   
   diff = lmean - rmean;
   return(diff);
}

NT_node insert_root(NT_node p, float diff)
{
   NT_node newp, prev, q, t;
   float dist, prevdist,td;

   
   newp = declare_tree_node( p->maxnseq);

   t = p->parent;
   prevdist = t->dist;

   p->parent = newp;

   dist = p->dist;

   p->dist = diff / 2;
   if (p->dist < 0.0) p->dist = 0.0;
   if (p->dist > dist) p->dist = dist;

   t->dist = dist - p->dist; 

   newp->left = t;
   newp->right = p;
   newp->parent = NULL;
   newp->dist = 0.0;
   newp->leaf = NODE;

   if (t->left == p) t->left = t->parent;
   else t->right = t->parent;

   prev = t;
   q = t->parent;

   t->parent = newp;

   while (q != NULL)
     {
        if (q->left == prev)
           {
              q->left = q->parent;
              q->parent = prev;
              td = q->dist;
              q->dist = prevdist;
              prevdist = td;
              prev = q;
              q = q->left;
           }
        else
           {
              q->right = q->parent;
              q->parent = prev;
              td = q->dist;
              q->dist = prevdist;
              prevdist = td;
              prev = q;
              q = q->right;
           }
    }

/*
   remove the old root node
*/
   q = prev;
   if (q->left == NULL)
      {
         dist = q->dist;
         q = q->right;
         q->dist += dist;
         q->parent = prev->parent;
         if (prev->parent->left == prev)
            prev->parent->left = q;
         else
            prev->parent->right = q;
         prev->right = NULL;
      }
   else
      {
         dist = q->dist;
         q = q->left;
         q->dist += dist;
         q->parent = prev->parent;
         if (prev->parent->left == prev)
            prev->parent->left = q;
         else
            prev->parent->right = q;
         prev->left = NULL;
      }

   return(newp);
}

/*********************************COPYRIGHT NOTICE**********************************/
/* Centre National de la Recherche Scientifique (CNRS) */
/*and */
/*Cedric Notredame */
/*Tue May 10 12:08:44     2005. */
/*All rights reserved.*/
/*This file is part of T-COFFEE.*/
/**/
/*    T-COFFEE is free software; you can redistribute it and/or modify*/
/*    it under the terms of the GNU General Public License as published by*/
/*    the Free Software Foundation; either version 2 of the License, or*/
/*    (at your option) any later version.*/
/**/
/*    T-COFFEE is distributed in the hope that it will be useful,*/
/*    but WITHOUT ANY WARRANTY; without even the implied warranty of*/
/*    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the*/
/*    GNU General Public License for more details.*/
/**/
/*    You should have received a copy of the GNU General Public License*/
/*    along with Foobar; if not, write to the Free Software*/
/*    Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA*/
/*...............................................                                                                                      |*/
/*  If you need some more information*/
/*  cedric.notredame@europe.com*/
/*...............................................                                                                                                                                     |*/
/**/
/**/
/*    */
/*********************************COPYRIGHT NOTICE**********************************/

Generated by  Doxygen 1.6.0   Back to index