~ubuntu-branches/ubuntu/utopic/paml/utopic

« back to all changes in this revision

Viewing changes to src/TreeTimeJeff.c

  • Committer: Bazaar Package Importer
  • Author(s): Pjotr Prins
  • Date: 2010-09-11 23:01:37 UTC
  • Revision ID: james.westby@ubuntu.com-20100911230137-jjf5d0blx5p0m9ba
Tags: upstream-4.4c
ImportĀ upstreamĀ versionĀ 4.4c

Show diffs side-by-side

added added

removed removed

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