-/*
- * treefunc.c
- *
- * (c) 2003 Gangolf Jobb (http://www.treefinder.de)
- *
- * Various data structures and methods for manipulating
- * and traversing trees
- * (e.g., to classify genes)
- *
- * This code may be distributed under the GNU GPL
- */
-
-#include <stdlib.h>
-#include <string.h>
-#include <math.h>
-
-/*============================================================================*/
-
-#define LABLEN 128
-
-typedef struct tree {
- char label[LABLEN];
- struct tree *branches,*next;
- double length;
- int mark;
-} Tree;
-
-/*....*/
-
-Tree *NewTree(void) {
- Tree *b;
-
- b=malloc(sizeof(Tree)); if(!b) return(NULL);
-
- *b->label='\0'; b->branches=NULL; b->next=NULL; b->length=0.0;
-
- return(b);
-}
-
-/*....*/
-
-void FreeTree(Tree *t) {
- Tree *b;
-
- while(t->branches) {b=t->branches; t->branches=b->next; FreeTree(b);}
-
- free(t);
-
-}
-
-/*============================================================================*/
-
-static Tree *current=NULL;
-
-static int tip_index;
-static int edge_index;
-static int node_index;
-
-enum {OK=0,ERROR}; /* one may invent more descriptive error codes */
-static int error=OK;
-
-/*============================================================================*/
-
-Tree *buildTreeFromPhylo_(
- int node,
- int *lowerNodes,
- int *upperNodes,
- double *edgeLengths,
- int nedges,
- char **tipLabels,
- int ntips
-) {
- Tree *t,*b,**bb;
- int i,j,n;
-
- t=NewTree();
-
- bb=&t->branches; n=0;
- for(i=0;i<nedges;i++) { if(lowerNodes[i]!=node) continue;
- j=upperNodes[i]; if(j==0) {error=ERROR; goto err;}
- if(j>0) { if(j>ntips) {error=ERROR; goto err;}
- b=NewTree(); strcpy(b->label,tipLabels[j-1]);
- } else { if(-j>nedges) {error=ERROR; goto err;}
- b=
- buildTreeFromPhylo_(j,lowerNodes,upperNodes,edgeLengths,nedges,tipLabels,ntips);
- }
- b->length=edgeLengths[i];
- *bb=b; bb=&b->next; n++;
- } if(n<2) {error=ERROR; goto err;}
- err:
- *bb=NULL;
-
- return(t);
-}
-
-/*....*/
-
-void buildTreeFromPhylo(
- int *lowerNodes,
- int *upperNodes,
- double *edgeLengths,
- int *nedges,
- char **tipLabels,
- int *ntips,
- int *result
-) {
-
- error=OK;
-
- if(current) {FreeTree(current); current=NULL;}
-
- if(*nedges<2||*ntips<2) {error=ERROR; *result=error; return;}
-
- current=buildTreeFromPhylo_(-1,lowerNodes,upperNodes,edgeLengths,*nedges,tipLabels,*ntips);
-
- if(error&¤t) {FreeTree(current); current=NULL;}
-
- *result=error;
-
-}
-
-/*============================================================================*/
-
-void destroyTree(int *result) {
-
- error=OK;
-
- if(current) {FreeTree(current); current=NULL;}
-
- *result=error;
-
-}
-
-/*============================================================================*/
-
-void getError(int *result) {
-
- *result=error;
-
- error=OK;
-
-}
-
-/*============================================================================*/
-
-int nTips_(Tree *t) {
- Tree *b;
- int n;
-
- if(!t->branches) return(1);
-
- n=0; for(b=t->branches;b;b=b->next) n+=nTips_(b);
-
- return(n);
-}
-
-/*....*/
-
-void nTips(int *result) {
-
- error=OK;
-
- if(!current) {error=ERROR; *result=0; return;}
-
- *result=nTips_(current);
-
-}
-
-/*============================================================================*/
-
-int nNodes_(Tree *t) {
- Tree *b;
- int n;
-
- if(!t->branches) return(1);
-
- n=1; for(b=t->branches;b;b=b->next) n+=nNodes_(b);
-
- return(n);
-}
-
-/*....*/
-
-void nNodes(int *result) {
-
- error=OK;
-
- if(!current) {error=ERROR; *result=0; return;}
-
- *result=nNodes_(current);
-
-}
-
-/*....*/
-
-void nEdges(int *result) {
-
- error=OK;
-
- if(!current) {error=ERROR; *result=0; return;}
-
- *result=nNodes_(current)-1;
-}
-
-/*============================================================================*/
-
-double markClasses_(Tree *t) {
- Tree *b; double destinct,sum;
-
- /* all tips above a marked ( == 1) node belong to the same class */
-
- if(!t->branches) {t->mark=1; return(t->length);}
-
- sum=0.; for(b=t->branches;b;b=b->next) sum+=markClasses_(b);
-
- destinct=nTips_(t)*(t->length); /* (t->length) == 0. at root */
-
- if(destinct>sum) { t->mark=1; return(destinct); }
-
- t->mark=0; return(sum);
-
-}
-
-/*....*/
-
-void getMisawaTajima__(Tree *t,int ignore,int *result) { /* maps tips to marked classes */
- Tree *b;
-
- if(t->mark&&!ignore) {node_index++; ignore=1;} /* marked nodes above a marked node will be ignored */
-
- if(!t->branches) {result[tip_index++]=node_index; return;}
-
- for(b=t->branches;b;b=b->next) getMisawaTajima__(b,ignore,result);
-
-}
-
-/*....*/
-
-void getMisawaTajima_(Tree *t,int *result) {
-
- markClasses_(t);
-
- tip_index=0; node_index=0;
-
- getMisawaTajima__(t,0,result);
-
-}
-
-/*....*/
-
-void getMisawaTajima(int *result) {
-
- error=OK;
-
- if(!current) {error=ERROR; return;}
-
- getMisawaTajima_(current,result);
-
-}