source: branches/stable/GDE/PHYLIP/penny.c

Last change on this file was 2175, checked in by westram, 20 years ago

upgrade to PHYLIP 3.6

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 22.7 KB
Line 
1
2#include "phylip.h"
3#include "disc.h"
4#include "wagner.h"
5
6/* version 3.6. (c) Copyright 1993-2002 by the University of Washington.
7   Written by Joseph Felsenstein, Akiko Fuseki, Sean Lamont, and Andrew Keeffe.
8   Permission is granted to copy and use this program provided no fee is
9   charged for it and provided that this copyright notice is not removed. */
10
11#define maxtrees        100  /* maximum number of trees to be printed out   */
12#define often           100  /* how often to notify how many trees examined */
13#define many            1000 /* how many multiples of howoften before stop  */
14
15typedef long *treenumbers;
16typedef double *valptr;
17typedef long *placeptr;
18
19#ifndef OLDC
20/* function prototypes */
21void   getoptions(void);
22void   allocrest(void);
23void   doinit(void);
24void   inputoptions(void);
25void   doinput(void);
26void   supplement(bitptr);
27void   evaluate(node2 *);
28void   addtraverse(node2 *,node2 *,node2 *,long *,long *,valptr,placeptr);
29void   addit(long);
30void   reroot(node2 *);
31void   describe(void);
32void   maketree(void);
33/* function prototypes */
34#endif
35
36Char infilename[FNMLNGTH], outfilename[FNMLNGTH], outtreename[FNMLNGTH], weightfilename[FNMLNGTH], ancfilename[FNMLNGTH], mixfilename[FNMLNGTH];
37node2 *root;
38long outgrno, rno, howmany, howoften, col, msets, ith;
39/*  outgrno indicates outgroup                                        */
40
41boolean weights, thresh, ancvar, questions, allsokal, allwagner,
42               mixture, simple, trout, noroot, didreroot, outgropt,
43               progress, treeprint, stepbox, ancseq, mulsets, firstset;
44boolean *ancone, *anczero, *ancone0, *anczero0, justwts;
45pointptr2 treenode;         /* pointers to all nodes in tree   */
46double fracdone, fracinc;
47double threshold;
48double *threshwt;
49bitptr wagner, wagner0;
50boolean *added;
51Char *guess;
52steptr numsteps;
53long **bestorders, **bestrees;
54steptr numsone, numszero;
55gbit *garbage;
56
57long examined, mults;
58boolean firsttime, done, full;
59double like, bestyet;
60treenumbers current, order;
61long fullset;
62bitptr zeroanc, oneanc;
63bitptr suppsteps;
64
65
66void getoptions()
67{
68  /* interactively set options */
69  long loopcount, loopcount2;
70  Char ch, ch2;
71
72  fprintf(outfile, "\nPenny algorithm, version %s\n",VERSION);
73  fprintf(outfile, " branch-and-bound to find all");
74  fprintf(outfile, " most parsimonious trees\n\n");
75  howoften = often;
76  howmany = many;
77  outgrno = 1;
78  outgropt = false;
79  simple = true;
80  thresh = false;
81  threshold = spp;
82  trout = true;
83  weights = false;
84  justwts = false;
85  ancvar = false;
86  allsokal = false;
87  allwagner = true;
88  mixture = false;
89  printdata = false;
90  progress = true;
91  treeprint = true;
92  stepbox = false;
93  ancseq = false;
94  loopcount = 0;
95  for(;;) {
96    cleerhome();
97    printf("\nPenny algorithm, version %s\n",VERSION);
98    printf(" branch-and-bound to find all most parsimonious trees\n\n");
99    printf("Settings for this run:\n");
100    printf("  X                     Use Mixed method?  %s\n",
101           mixture ? "Yes" : "No");
102    printf("  P                     Parsimony method?  %s\n",
103           (allwagner && !mixture)  ? "Wagner"       :
104           (!(allwagner || mixture)) ? "Camin-Sokal"  : "(methods in mixture)");
105    printf("  F        How often to report, in trees:%5ld\n",howoften);
106    printf("  H        How many groups of%5ld trees:%6ld\n",howoften,howmany);
107    printf("  O                        Outgroup root?");
108    if (outgropt)
109      printf("  Yes, at species number%3ld\n", outgrno);
110    else
111      printf("  No, use as outgroup species%3ld\n", outgrno);
112    printf("  S           Branch and bound is simple?  %s\n",
113           simple ? "Yes" : "No. reconsiders order of species");
114    printf("  T              Use Threshold parsimony?");
115    if (thresh)
116      printf("  Yes, count steps up to%4.1f per char.\n", threshold);
117    else
118      printf("  No, use ordinary parsimony\n");
119    printf("  A   Use ancestral states in input file?  %s\n",
120           ancvar ? "Yes" : "No");
121    printf("  W                       Sites weighted?  %s\n",
122           (weights ? "Yes" : "No"));
123    printf("  M           Analyze multiple data sets?");
124    if (mulsets)
125        printf("  Yes, %2ld %s\n", msets,
126               (justwts ? "sets of weights" : "data sets"));
127    else
128      printf("  No\n");
129    printf("  0   Terminal type (IBM PC, ANSI, none)?  %s\n",
130           ibmpc ? "IBM PC" : ansi  ? "ANSI" : "(none)");
131    printf("  1    Print out the data at start of run  %s\n",
132           printdata ? "Yes" : "No");
133    printf("  2  Print indications of progress of run  %s\n",
134           progress ? "Yes" : "No");
135    printf("  3                        Print out tree  %s\n",
136           treeprint ? "Yes" : "No");
137    printf("  4     Print out steps in each character  %s\n",
138           stepbox ? "Yes" : "No");
139    printf("  5     Print states at all nodes of tree  %s\n",
140           ancseq ? "Yes" : "No");
141    printf("  6       Write out trees onto tree file?  %s\n",
142           trout ? "Yes" : "No");
143    if(weights && justwts){
144        printf(
145         "WARNING:  W option and Multiple Weights options are both on.  ");
146        printf(
147         "The W menu option is unnecessary and has no additional effect. \n");
148    }
149    printf("\nAre these settings correct?");
150    printf(" (type Y or the letter for one to change)\n");
151#ifdef WIN32
152    phyFillScreenColor();
153#endif
154    scanf("%c%*[^\n]", &ch);
155    getchar();
156    uppercase(&ch);
157    if (ch == 'Y')
158      break;
159    if (strchr("WHFSOMPATX1234560",ch) != NULL){
160      switch (ch) {
161       
162      case 'X':
163        mixture = !mixture;
164        break;
165       
166      case 'P':
167        allwagner = !allwagner;
168        break;
169       
170      case 'A':
171        ancvar = !ancvar;
172        break;
173       
174      case 'H':
175        inithowmany(&howmany, howoften);
176        break;
177       
178      case 'F':
179        inithowoften(&howoften);
180        break;
181       
182      case 'S':
183        simple = !simple;
184        break;
185       
186      case 'O':
187        outgropt = !outgropt;
188        if (outgropt)
189          initoutgroup(&outgrno, spp);
190        else outgrno = 1;
191        break;
192       
193      case 'T':
194        thresh = !thresh;
195        if (thresh)
196          initthreshold(&threshold);
197        break;
198
199      case 'W':
200        weights = !weights;
201        break;
202
203      case 'M':
204        mulsets = !mulsets;
205        if (mulsets){
206            printf("Multiple data sets or multiple weights?");
207          loopcount2 = 0;
208          do {
209            printf(" (type D or W)\n");
210#ifdef WIN32
211            phyFillScreenColor();
212#endif
213            scanf("%c%*[^\n]", &ch2);
214            getchar();
215            if (ch2 == '\n')
216              ch2 = ' ';
217            uppercase(&ch2);
218            countup(&loopcount2, 10);
219          } while ((ch2 != 'W') && (ch2 != 'D'));
220          justwts = (ch2 == 'W');
221          if (justwts)
222            justweights(&msets);
223          else
224            initdatasets(&msets);
225        }
226        break;
227       
228      case '0':
229        initterminal(&ibmpc, &ansi);
230        break;
231       
232      case '1':
233        printdata = !printdata;
234        break;
235       
236      case '2':
237        progress = !progress;
238        break;
239       
240      case '3':
241        treeprint = !treeprint;
242        break;
243       
244      case '4':
245        stepbox = !stepbox;
246        break;
247       
248      case '5':
249        ancseq = !ancseq;
250        break;
251       
252      case '6':
253        trout = !trout;
254        break;
255      }
256    } else
257      printf("Not a possible option!\n");
258    countup(&loopcount, 100);
259  }
260  allsokal = (!allwagner && !mixture);
261}  /* getoptions */
262
263
264void allocrest()
265{
266  long i;
267
268  weight = (steptr)Malloc(chars*sizeof(steptr));
269  threshwt = (double *)Malloc(chars*sizeof(double));
270  bestorders = (long **)Malloc(maxtrees*sizeof(long *));
271  bestrees = (long **)Malloc(maxtrees*sizeof(long *));
272  for (i = 1; i <= maxtrees; i++) {
273    bestorders[i - 1] = (long *)Malloc(spp*sizeof(long));
274    bestrees[i - 1] = (long *)Malloc(spp*sizeof(long));
275  }
276  numsteps = (steptr)Malloc(chars*sizeof(steptr));
277  guess = (Char *)Malloc(chars*sizeof(Char));
278  numszero = (steptr)Malloc(chars*sizeof(steptr));
279  numsone = (steptr)Malloc(chars*sizeof(steptr));
280  current = (treenumbers)Malloc(spp*sizeof(long));
281  order = (treenumbers)Malloc(spp*sizeof(long));
282  nayme = (naym *)Malloc(spp*sizeof(naym));
283  added = (boolean *)Malloc(nonodes*sizeof(boolean));
284  ancone = (boolean *)Malloc(chars*sizeof(boolean));
285  anczero = (boolean *)Malloc(chars*sizeof(boolean));
286  ancone0 = (boolean *)Malloc(chars*sizeof(boolean));
287  anczero0 = (boolean *)Malloc(chars*sizeof(boolean));
288  wagner = (bitptr)Malloc(words*sizeof(long));
289  wagner0 = (bitptr)Malloc(words*sizeof(long));
290  zeroanc = (bitptr)Malloc(words*sizeof(long));
291  oneanc = (bitptr)Malloc(words*sizeof(long));
292  suppsteps = (bitptr)Malloc(words*sizeof(long));
293  extras = (steptr)Malloc(chars*sizeof(steptr));
294}  /* allocrest */
295
296
297void doinit()
298{
299  /* initializes variables */
300
301  inputnumbers(&spp, &chars, &nonodes, 1);
302  words = chars / bits + 1;
303  getoptions();
304  if (printdata)
305    fprintf(outfile, "%2ld species, %3ld characters\n", spp, chars);
306  alloctree2(&treenode);
307  setuptree2(treenode);
308  allocrest();
309}  /* doinit */
310
311void inputoptions()
312{
313  /* input the information on the options */
314  long i;
315  if(justwts){
316      if(firstset){
317          scan_eoln(infile);
318          if (ancvar) {
319              inputancestorsnew(anczero0, ancone0);
320          }
321          if (mixture) {
322              inputmixturenew(wagner0);
323          }
324      }
325      for (i = 0; i < (chars); i++)
326          weight[i] = 1;
327      inputweights(chars, weight, &weights);
328      for (i = 0; i < (words); i++) {
329          if (mixture)
330              wagner[i] = wagner0[i];
331          else if (allsokal)
332              wagner[i] = 0;
333          else
334              wagner[i] = (1L << (bits + 1)) - (1L << 1);
335      }
336  }
337  else {
338      if (!firstset) {
339          samenumsp(&chars, ith);
340      }
341      scan_eoln(infile);
342      for (i = 0; i < (chars); i++)
343          weight[i] = 1;
344      if (ancvar) {
345          inputancestorsnew(anczero0, ancone0);
346      }
347      if (mixture) {
348          inputmixturenew(wagner0);
349      }
350      if (weights)
351          inputweights(chars, weight, &weights);
352      for (i = 0; i < (words); i++) {
353          if (mixture)
354              wagner[i] = wagner0[i];
355          else if (allsokal)
356              wagner[i] = 0;
357          else
358              wagner[i] = (1L << (bits + 1)) - (1L << 1);
359      }
360  }
361  for (i = 0; i < (chars); i++) {
362      if (!ancvar) {
363          anczero[i] = true;
364          ancone[i] = (((1L << (i % bits + 1)) & wagner[i / bits]) != 0);
365      } else {
366          anczero[i] = anczero0[i];
367          ancone[i] = ancone0[i];
368      }
369  }
370  noroot = true;
371  questions = false;
372  for (i = 0; i < (chars); i++) {
373      if (weight[i] > 0) {
374          noroot = (noroot && ancone[i] && anczero[i] &&
375                    ((((1L << (i % bits + 1)) & wagner[i / bits]) != 0)
376                     || threshold <= 2.0));
377      }
378      questions = (questions || (ancone[i] && anczero[i]));
379      threshwt[i] = threshold * weight[i];
380  }
381}  /* inputoptions */
382
383
384void doinput()
385{
386  /* reads the input data */
387  inputoptions();
388  if(!justwts || firstset)
389      inputdata2(treenode);
390}  /* doinput */
391
392
393void supplement(bitptr suppsteps)
394{
395  /* determine minimum number of steps more which will
396     be added when rest of species are put in tree */
397  long i, j, k, l;
398  long defone, defzero, a;
399
400  k = 0;
401  for (i = 0; i < (words); i++) {
402    defone = 0;
403    defzero = 0;
404    a = 0;
405    for (l = 1; l <= bits; l++) {
406      k++;
407      if (k <= chars) {
408        if (!ancone[k - 1])
409          defzero = ((long)defzero) | (1L << l);
410        if (!anczero[k - 1])
411          defone = ((long)defone) | (1L << l);
412      }
413    }
414    for (j = 0; j < (spp); j++) {
415      defone |= treenode[j]->empstte1[i] & (~treenode[j]->empstte0[i]);
416      defzero |= treenode[j]->empstte0[i] & (~treenode[j]->empstte1[i]);
417      if (added[j])
418        a |= defone & defzero;
419    }
420    suppsteps[i] = defone & defzero & (~a);
421  }
422}  /* supplement */
423
424
425void evaluate(node2 *r)
426{
427  /* Determines the number of steps needed for a tree.
428     This is the minimum number needed to evolve chars on
429     this tree */
430  long i, stepnum, smaller;
431  double sum;
432
433  sum = 0.0;
434  for (i = 0; i < (chars); i++) {
435    numszero[i] = 0;
436    numsone[i] = 0;
437  }
438  supplement(suppsteps);
439  for (i = 0; i < (words); i++)
440    zeroanc[i] =fullset;
441  full = true;
442  postorder(r, fullset, full, wagner, zeroanc);
443  cpostorder(r, full, zeroanc, numszero, numsone);
444  count(r->fulstte1, zeroanc, numszero, numsone);
445  count(suppsteps, zeroanc, numszero, numsone);
446  for (i = 0; i < (words); i++)
447    zeroanc[i] = 0;
448  full = false;
449  postorder(r, fullset, full, wagner, zeroanc);
450  cpostorder(r, full, zeroanc, numszero, numsone);
451  count(r->empstte0, zeroanc, numszero, numsone);
452  count(suppsteps, zeroanc, numszero, numsone);
453  for (i = 0; i < (chars); i++) {
454    smaller = spp * weight[i];
455    numsteps[i] = smaller;
456    if (anczero[i]) {
457      numsteps[i] = numszero[i];
458      smaller = numszero[i];
459    }
460    if (ancone[i] && numsone[i] < smaller)
461      numsteps[i] = numsone[i];
462    stepnum = numsteps[i] + extras[i];
463    if (stepnum <= threshwt[i])
464      sum += stepnum;
465    else
466      sum += threshwt[i];
467    guess[i] = '?';
468    if (!ancone[i] || (anczero[i] && numszero[i] < numsone[i]))
469      guess[i] = '0';
470    else if (!anczero[i] || (ancone[i] && numsone[i] < numszero[i]))
471      guess[i] = '1';
472  }
473  if (examined == 0 && mults == 0)
474    bestyet = -1.0;
475  like = sum;
476}  /* evaluate */
477
478
479void addtraverse(node2 *a, node2 *b, node2 *c, long *m, long *n,
480                        valptr valyew, placeptr place)
481{
482  /* traverse all places to add b */
483  if (done)
484    return;
485  if ((*m) <= 2 || !(noroot && (a == root || a == root->next->back))) {
486    add3(a, b, c, &root, treenode);
487    (*n)++;
488    evaluate(root);
489    examined++;
490    if (examined == howoften) {
491      examined = 0;
492      mults++;
493      if (mults == howmany)
494        done = true;
495      if (progress) {
496        printf("%6ld", mults);
497        if (bestyet >= 0)
498          printf("%18.5f", bestyet);
499        else
500          printf("         -        ");
501        printf("%17ld%20.2f\n", nextree - 1, fracdone * 100);
502#ifdef WIN32
503        phyFillScreenColor();
504#endif
505      }
506    }
507    valyew[(*n) - 1] = like;
508    place[(*n) - 1] = a->index;
509    re_move3(&b, &c, &root, treenode);
510  }
511  if (!a->tip) {
512    addtraverse(a->next->back, b, c, m,n,valyew,place);
513    addtraverse(a->next->next->back, b, c, m,n,valyew,place);
514  }
515}  /* addtraverse */
516
517
518void addit(long m)
519{
520  /* adds the species one by one, recursively */
521  long n;
522  valptr valyew;
523  placeptr place;
524  long i, j, n1, besttoadd = 0;
525  valptr bestval;
526  placeptr bestplace;
527  double oldfrac, oldfdone, sum, bestsum;
528
529  valyew = (valptr)Malloc(nonodes*sizeof(double));
530  bestval = (valptr)Malloc(nonodes*sizeof(double));
531  place = (placeptr)Malloc(nonodes*sizeof(long));
532  bestplace = (placeptr)Malloc(nonodes*sizeof(long));
533  if (simple && !firsttime) {
534    n = 0;
535    added[order[m - 1] - 1] = true;
536    addtraverse(root, treenode[order[m - 1] - 1],
537                treenode[spp + m - 2], &m,&n,valyew,place);
538    besttoadd = order[m - 1];
539    memcpy(bestplace, place, nonodes*sizeof(long));
540    memcpy(bestval, valyew, nonodes*sizeof(double));
541  } else {
542    bestsum = -1.0;
543    for (i = 1; i <= (spp); i++) {
544      if (!added[i - 1]) {
545        n = 0;
546        added[i - 1] = true;
547        addtraverse(root, treenode[i - 1], treenode[spp + m - 2], &m,
548                    &n,valyew,place);
549        added[i - 1] = false;
550        sum = 0.0;
551        for (j = 0; j < (n); j++)
552          sum += valyew[j];
553        if (sum > bestsum) {
554          bestsum = sum;
555          besttoadd = i;
556          memcpy(bestplace, place, nonodes*sizeof(long));
557          memcpy(bestval, valyew, nonodes*sizeof(double));
558        }
559      }
560    }
561  }
562  order[m - 1] = besttoadd;
563  memcpy(place, bestplace, nonodes*sizeof(long));
564  memcpy(valyew, bestval, nonodes*sizeof(double));
565  shellsort(valyew, place, n);
566  oldfrac = fracinc;
567  oldfdone = fracdone;
568  n1 = 0;
569  for (i = 0; i < (n); i++) {
570    if (valyew[i] <= bestyet || bestyet < 0.0)
571      n1++;
572  }
573  if (n1 > 0)
574    fracinc /= n1;
575  for (i = 0; i < (n); i++) {
576    if (valyew[i] <= bestyet || bestyet < 0.0) {
577      current[m - 1] = place[i];
578      add3(treenode[place[i] - 1], treenode[besttoadd - 1],
579          treenode[spp + m - 2], &root, treenode);
580      added[besttoadd - 1] = true;
581      if (m < spp)
582        addit(m + 1);
583      else {
584        if (valyew[i] < bestyet || bestyet < 0.0) {
585          nextree = 1;
586          bestyet = valyew[i];
587        }
588        if (nextree <= maxtrees) {
589          memcpy(bestorders[nextree - 1], order,
590                 spp*sizeof(long));
591          memcpy(bestrees[nextree - 1], current,
592                 spp*sizeof(long));
593        }
594        nextree++;
595        firsttime = false;
596      }
597      re_move3(&treenode[besttoadd - 1], &treenode[spp + m - 2], &root, treenode);
598      added[besttoadd - 1] = false;
599    }
600    fracdone += fracinc;
601  }
602  fracinc = oldfrac;
603  fracdone = oldfdone;
604  free(valyew);
605  free(bestval);
606  free(place);
607  free(bestplace);
608}  /* addit */
609
610
611void reroot(node2 *outgroup)
612{
613  /* reorients tree, putting outgroup in desired position. */
614  node2 *p, *q, *newbottom, *oldbottom;
615
616  if (outgroup->back->index == root->index)
617    return;
618  newbottom = outgroup->back;
619  p = treenode[newbottom->index - 1]->back;
620  while (p->index != root->index) {
621    oldbottom = treenode[p->index - 1];
622    treenode[p->index - 1] = p;
623    p = oldbottom->back;
624  }
625  p = root->next;
626  q = root->next->next;
627  p->back->back = q->back;
628  q->back->back = p->back;
629  p->back = outgroup;
630  q->back = outgroup->back;
631  outgroup->back->back = root->next->next;
632  outgroup->back = root->next;
633  treenode[newbottom->index - 1] = newbottom;
634}  /* reroot */
635
636
637void describe()
638{
639  /* prints ancestors, steps and table of numbers of steps in
640     each character */
641
642  if (stepbox) {
643    putc('\n', outfile);
644    writesteps(weights, numsteps);
645  }
646  if (questions && (!noroot || didreroot))
647    guesstates(guess);
648  if (ancseq) {
649    hypstates(fullset, full, noroot, didreroot, root, wagner,
650                zeroanc, oneanc, treenode, guess, garbage);
651    putc('\n', outfile);
652  }
653  if (trout) {
654    col = 0;
655    treeout2(root, &col, root);
656  }
657}  /* describe */
658
659
660void maketree()
661{
662  /* tree construction recursively by branch and bound */
663  long i, j, k;
664  node2 *dummy;
665
666  fullset = (1L << (bits + 1)) - (1L << 1);
667  if (progress) {
668    printf("\nHow many\n");
669    printf("trees looked                                       Approximate\n");
670    printf("at so far      Length of        How many           percentage\n");
671    printf("(multiples     shortest tree    trees this long    searched\n");
672    printf("of %4ld):      found so far     found so far       so far\n",
673           howoften);
674    printf("----------     ------------     ------------       ------------\n");
675#ifdef WIN32
676    phyFillScreenColor();
677#endif
678  }
679  done = false;
680  mults = 0;
681  examined = 0;
682  nextree = 1;
683  root = treenode[0];
684  firsttime = true;
685  for (i = 0; i < (spp); i++)
686    added[i] = false;
687  added[0] = true;
688  order[0] = 1;
689  k = 2;
690  fracdone = 0.0;
691  fracinc = 1.0;
692  bestyet = -1.0;
693  addit(k);
694  if (done) {
695    if (progress) {
696      printf("Search broken off!  Not guaranteed to\n");
697      printf(" have found the most parsimonious trees.\n");
698    }
699    if (treeprint) {
700      fprintf(outfile, "Search broken off!  Not guaranteed to\n");
701      fprintf(outfile, " have found the most parsimonious\n");
702      fprintf(outfile, " trees, but here is what we found:\n");
703    }
704  }
705  if (treeprint) {
706    fprintf(outfile, "\nrequires a total of %18.3f\n\n", bestyet);
707    if (nextree == 2)
708      fprintf(outfile, "One most parsimonious tree found:\n");
709    else
710      fprintf(outfile, "%5ld trees in all found\n", nextree - 1);
711  }
712  if (nextree > maxtrees + 1) {
713    if (treeprint)
714      fprintf(outfile, "here are the first%4ld of them\n", (long)maxtrees);
715    nextree = maxtrees + 1;
716  }
717  if (treeprint)
718    putc('\n', outfile);
719  for (i = 0; i < (spp); i++)
720    added[i] = true;
721  for (i = 0; i <= (nextree - 2); i++) {
722    for (j = k; j <= (spp); j++)
723      add3(treenode[bestrees[i][j - 1] - 1],
724          treenode[bestorders[i][j - 1] - 1], treenode[spp + j - 2],
725          &root, treenode);
726    if (noroot)
727      reroot(treenode[outgrno - 1]);
728    didreroot = (outgropt && noroot);
729    evaluate(root);
730    printree(treeprint, noroot, didreroot, root);
731    describe();
732    for (j = k - 1; j < (spp); j++)
733      re_move3(&treenode[bestorders[i][j] - 1], &dummy, &root, treenode);
734  }
735  if (progress) {
736    printf("\nOutput written to file \"%s\"\n\n", outfilename);
737    if (trout)
738      printf("Trees also written onto file \"%s\"\n\n", outtreename);
739  }
740  if (ancseq)
741    freegarbage(&garbage);
742}  /* maketree */
743
744
745int main(int argc, Char *argv[])
746{  /* Penny's branch-and-bound method */
747   /* Reads in the number of species, number of characters,
748     options and data.  Then finds all most parsimonious trees */
749#ifdef MAC
750  argc = 1;                /* macsetup("Penny","");                */
751  argv[0] = "Penny";
752#endif
753  init(argc,argv);
754  openfile(&infile,INFILE,"input file", "r",argv[0],infilename);
755  openfile(&outfile,OUTFILE,"output file", "w",argv[0],outfilename);
756  ibmpc = IBMCRT;
757  ansi = ANSICRT;
758  mulsets = false;
759  msets = 1;
760  firstset = true;
761  garbage = NULL;
762  bits = 8*sizeof(long) - 1;
763  doinit();
764  if (weights || justwts)
765    openfile(&weightfile,WEIGHTFILE,"weights file","r",argv[0],weightfilename);
766  if (trout)
767    openfile(&outtree,OUTTREE,"output tree file", "w",argv[0],outtreename);
768  if(ancvar)
769      openfile(&ancfile,ANCFILE,"ancestors file", "r",argv[0],ancfilename);
770  if(mixture)
771      openfile(&mixfile,MIXFILE,"mixture file", "r",argv[0],mixfilename);
772 
773  for (ith = 1; ith <= msets; ith++) {
774    if(firstset){
775        if (allsokal && !mixture)
776            fprintf(outfile, "Camin-Sokal parsimony method\n\n");
777        if (allwagner && !mixture)
778            fprintf(outfile, "Wagner parsimony method\n\n");
779    }
780    doinput();
781    if (msets > 1 && !justwts) {
782      fprintf(outfile, "Data set # %ld:\n\n",ith);
783      if (progress)
784        printf("\nData set # %ld:\n",ith);
785    } 
786    if (justwts){
787        if(firstset && mixture && printdata)
788            printmixture(outfile, wagner);
789        fprintf(outfile, "Weights set # %ld:\n\n", ith);
790        if (progress)
791            printf("\nWeights set # %ld:\n\n", ith);
792    }
793    else if (mixture && printdata)
794        printmixture(outfile, wagner);
795    if (printdata){
796        if (weights || justwts)
797            printweights(outfile, 0, chars, weight, "Characters");
798        if (ancvar)
799            printancestors(outfile, anczero, ancone);
800    }
801    if (ith == 1)
802      firstset = false;
803    maketree();
804  }
805  FClose(infile);
806  FClose(outfile);
807  FClose(outtree);
808#ifdef MAC
809  fixmacfile(outfilename);
810  fixmacfile(outtreename);
811#endif
812#ifdef WIN32
813  phyRestoreConsoleAttributes();
814#endif
815  return 0;
816}  /* Penny's branch-and-bound method */
817
Note: See TracBrowser for help on using the repository browser.