src/TreeTimeJeff.c

   1 /* TreeTimeJeff.c\r
   2    Copyright, Ziheng Yang, June 2003.\r
   3 \r
   4      cc -O3 -o TreeTimeJeff TreeTimeJeff.c tools.c -lm\r
   5      cl -Ot -O2 TreeTimeJeff.c tools.c\r
   6 \r
   7      TreeTimeJeff <MultidivtimeOutputFile>\r
   8 \r
   9    This reads the output from Jeff Thorne's multidivtime to print out\r
  10    a tree with branch lengthes calculated using the posterior means of\r
  11    divergence times.  The output tree can then be read in from Rod\r
  12    Page's TreeView.  This is for making a tree like Figure 1 of Yang &\r
  13    Yoder (2003 Syst Biol), the cute-looking paper on mouse lemurs.\r
  14    Suppose the output from multidivtime is o.md.  Run the program by\r
  15 \r
  16    TreeTimeJeff o.md\r
  17    TreeTimeJeff o.md > z.trees\r
  18 \r
  19    Then edit out everything except the tree in z.trees.  Read the tree\r
  20    from TreeView, change fonts etc., save it as a picture file and\r
  21    import it into powerpoint.  Use excel to make up an x-axis, of\r
  22    reasonable length, to be used as the time scale.  Delete everything\r
  23    else and copy the x-axis into powerpoint.  Resize the tree to match\r
  24    up with the time axis.\r
  25 \r
  26 */\r
  27 \r
  28 #include "paml.h"\r
  29 #define NS            1000\r
  30 #define NBRANCH       (NS*2-2)\r
  31 #define MAXNSONS      100\r
  32 \r
  33 struct CommonInfo {\r
  34    char *z[2*NS-1], spname[NS][20], cleandata;\r
  35    int ns, ls, npatt, np, ntime, ncode, clock, model, icode;\r
  36    int seqtype, *pose;\r
  37    double *conP, *fpatt;\r
  38 }  com;\r
  39 struct TREEB {\r
  40    int nbranch, nnode, root, branches[NBRANCH][2];\r
  41 }  tree;\r
  42 struct TREEN {\r
  43    int father, nson, sons[NS], ibranch;\r
  44    double branch, age, label, *conP;\r
  45   char *nodeStr;\r
  46 }  nodes[2*NS-1];\r
  47 \r
  48 #define NODESTRUCTURE\r
  49 #include "treesub.c"\r
  50 \r
  51 void main (int argc, char*argv[])\r
  52 {\r
  53    int  lline=32000, i,j, ch, jeffnode, inode;\r
  54    char line[32000], mcmcf[96]="o.multidivtime";\r
  55    FILE *fmcmc;\r
  56    double t;\r
  57 \r
  58    puts("Usage:\n\tTreeTimeJeff <MultidivtimeOutputFile>\n");\r
  59    if(argc>1) strcpy(mcmcf, argv[1]);\r
  60    fmcmc=gfopen(mcmcf, "r");\r
  61 \r
  62    /* Read root node number */\r
  63    for( ; ; ) {\r
  64            if(fgets(line, lline, fmcmc) == NULL) error2("EOF mcmc file");\r
  65       if(strstr(line, "Root node number of master tree is")==NULL) continue;\r
  66       sscanf(line+35, "%d", &j);\r
  67       com.ns=j/2+1;\r
  68       break;\r
  69    }\r
  70    printf("Tree has %d taxa.\n\n", com.ns);\r
  71 \r
  72    /* read tree.  JeffNode read into [].branch */\r
  73    for(; ; ) {\r
  74            ch=fgetc(fmcmc);\r
  75            if(ch==EOF) error2("EOF treefile");\r
  76            if(ch=='(') \r
  77          { ungetc(ch,fmcmc); break; }\r
  78    }\r
  79    ReadTreeN(fmcmc, &i, &j, 2, 0);\r
  80    OutTreeN(F0,1,0);  FPN(F0);  FPN(F0);\r
  81 \r
  82    /* read posterior time estimates */\r
  83    for(i=0; i<tree.nnode; i++) nodes[i].age=0;\r
  84    for( ; ; ) {\r
  85            if(fgets(line, lline, fmcmc) == NULL) error2("EOF mcmc file");\r
  86       if(strstr(line, "Actual time node")==NULL) continue;\r
  87       sscanf(line+17, "%d =%lf", &jeffnode, &t);\r
  88       if(jeffnode<com.ns) {\r
  89          if(t>0) nodes[inode].age=t;\r
  90       }\r
  91       else {\r
  92          inode=tree.nnode-1+com.ns-jeffnode;\r
  93          printf("JeffNode %3d ZihengNode %3d time %9.6f\n", jeffnode,inode+1,t);\r
  94          nodes[inode].age=t;\r
  95          if(inode==com.ns) break;\r
  96          if(jeffnode-nodes[inode].branch != 0)\r
  97             printf(" node number error. ");\r
  98       }\r
  99    }\r
 100    for(i=0; i<tree.nnode; i++) \r
 101       if(i!=tree.root) nodes[i].branch=nodes[nodes[i].father].age-nodes[i].age;\r
 102 \r
 103    FPN(F0);  OutTreeN(F0,1,1);  FPN(F0);\r
 104    fclose(fmcmc);\r
 105    exit(0);\r
 106 }\r