newickstring.cpp
上传用户:yhdzpy8989
上传日期:2007-06-13
资源大小:13604k
文件大小:11k
源码类别:

生物技术

开发平台:

C/C++

  1. /*
  2.  * ===========================================================================
  3.  * PRODUCTION $Log: newickstring.cpp,v $
  4.  * PRODUCTION Revision 1000.1  2004/06/01 18:10:14  gouriano
  5.  * PRODUCTION PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.2
  6.  * PRODUCTION
  7.  * ===========================================================================
  8.  */
  9. /*  $Id: newickstring.cpp,v 1000.1 2004/06/01 18:10:14 gouriano Exp $
  10. * ===========================================================================
  11. *
  12. *                            PUBLIC DOMAIN NOTICE
  13. *               National Center for Biotechnology Information
  14. *
  15. *  This software/database is a "United States Government Work" under the
  16. *  terms of the United States Copyright Act.  It was written as part of
  17. *  the author's official duties as a United States Government employee and
  18. *  thus cannot be copyrighted.  This software/database is freely available
  19. *  to the public for use. The National Library of Medicine and the U.S.
  20. *  Government have not placed any restriction on its use or reproduction.
  21. *
  22. *  Although all reasonable efforts have been taken to ensure the accuracy
  23. *  and reliability of the software and data, the NLM and the U.S.
  24. *  Government do not and cannot warrant the performance or results that
  25. *  may be obtained by using this software or data. The NLM and the U.S.
  26. *  Government disclaim all warranties, express or implied, including
  27. *  warranties of performance, merchantability or fitness for any particular
  28. *  purpose.
  29. *
  30. *  Please cite the author in any work or product based on this material.
  31. *
  32. * ===========================================================================
  33. *
  34. * Author:  Richard Desper
  35. *
  36. * File Description:  newickstring.cpp
  37. *
  38. *    A part of the Miminum Evolution algorithm
  39. *
  40. */
  41. #include <ncbi_pch.hpp>
  42. #include <stdio.h>
  43. #include <stdlib.h>
  44. #include <math.h>
  45. #include "graph.h"
  46. #include "fastme.h"
  47. #include "newick.h"
  48. /*
  49. Hi Chris,
  50. I've written a function which will write the meTree as a string.  It's in
  51. the newickstring.c file, called NewickPrintTreeString.
  52. I needed to add one constant to the newick.h file to do this, MAX_DIGITS.
  53. I'm assuming 20 digits would suffice - if you need more you'll need to change 
  54. this value, or the program may not allot a sufficiently long string.  
  55. Rick
  56. */
  57. BEGIN_NCBI_SCOPE
  58. BEGIN_SCOPE(fastme)
  59. boolean leaf(meNode *v);
  60. boolean whiteSpace(char c)
  61. {
  62.   if ((' ' == c) || ('t' == c) || ('n' == c))
  63.     return(TRUE_FASTME);
  64.   else return(FALSE_FASTME);
  65. }
  66. /*decodeNewickSubmeTree is used to turn a string of the form
  67.   "(v1:d1,v2:d2,(subtree) v3:d3....vk:dk) subroot:d," into a subtree
  68.   rooted at subrooted, with corresponding subtrees and leaves at v1
  69.   through vk.  It is called recursively on subtrees*/
  70. meNode *decodeNewickSubtree(char *treeString, meTree *T, int *uCount)
  71. {
  72.   meNode *thisNode = NULL;
  73.   meNode *centerNode;
  74.   double thisWeight;
  75.   meEdge *thisEdge;
  76.   char label[MAX_LABEL_LENGTH];
  77.   char stringWeight[MAX_LABEL_LENGTH];
  78.   int state;
  79.   int i = 0;
  80.   int j;
  81.   int left,right;
  82.   int parcount;
  83.   sprintf(label,"Default Label");
  84.   left = right = 0;
  85.   parcount = 0;
  86.   state = ReadOpenParenthesis;
  87.   if('(' == treeString[0])
  88.     parcount++;
  89.   centerNode = makeNode(label,NULL,nodeCount++);
  90.   T->size++;
  91.   while(parcount > 0)
  92.     {
  93.       while(whiteSpace(treeString[i]))
  94. i++;
  95.       switch(state) 
  96. {
  97. case(ReadOpenParenthesis):
  98.   if('(' != treeString[0])
  99.     {
  100.       fprintf(stderr,"Error reading subtree.n");
  101.       exit(EXIT_FAILURE);
  102.     }
  103.   i++;
  104.   state = ReadSubTree;
  105.   break;
  106. case(ReadSubTree):
  107.   if('(' == treeString[i])  /*if treeString[i] is a left parenthesis,
  108.       we scan down the string until we find its partner.
  109.       the relative number of '('s vs. ')'s is counted
  110.       by the variable parcount*/
  111.     {
  112.       left = i++;
  113.       parcount++;
  114.       while(parcount > 1)
  115. {
  116.   while (('(' != treeString[i]) && (')' != treeString[i]))
  117.     i++;  /*skip over characters which are not parentheses*/
  118.   if('(' == treeString[i])
  119.     parcount++;
  120.   else if (')' == treeString[i])
  121.     parcount--;
  122.   i++;
  123. }  /*end while */
  124.       right = i;  /*at this point, the submeTree string goes from 
  125.     treeString[left] to treeString[right - 1]*/
  126.       thisNode = decodeNewickSubtree(treeString + left,T,uCount);  /*note that this
  127.       step will put 
  128.       thisNode in T*/
  129.       i = right;  /*having created the meNode for the subtree, we move
  130.     to assigning the label for the new node.
  131.     treeString[right] will be the start of this label */
  132.     } /* end if ('(' == treeString[i]) condition */
  133.   else
  134.     {   
  135.       thisNode = makeNode(label,NULL,nodeCount++);
  136.       T->size++;
  137.     }
  138.   state = ReadLabel;
  139.   break;
  140. case(ReadLabel):
  141.   left = i;  /*recall "left" is the left marker for the substring, "right" the right*/
  142.   if (':' == treeString[i])   /*this means an internal node?*/
  143.     {
  144.       sprintf(thisNode->label,"I%d",(*uCount)++);
  145.       right = i;
  146.     }
  147.   else
  148.     {
  149.       while((':' != treeString[i]) && (',' != treeString[i]))
  150. i++;
  151.       right = i;
  152.       j = 0;
  153.       for(i = left; i < right;i++)
  154. if(!(whiteSpace(treeString[i])))
  155.   thisNode->label[j++] = treeString[i];
  156.       thisNode->label[j] = '';
  157.     }       
  158.   if(':' == treeString[right])
  159.     state = ReadWeight;
  160.   else
  161.     state = ReadSubTree;
  162.   i = right + 1;
  163.   break;
  164. case(ReadWeight):
  165.   left = i;
  166.   while
  167.     ((',' != treeString[i]) && (')' != treeString[i]))
  168.     i++;
  169.   right = i;
  170.   if (',' == treeString[right])
  171.     state = ReadSubTree;
  172.   else
  173.     parcount--;
  174.   j = 0;
  175.   for(i = left; i < right; i++)
  176.     stringWeight[j++] = treeString[i];
  177.   stringWeight[j] = '';
  178.   thisWeight = atof(stringWeight);
  179.   thisEdge = makeEdge(label,centerNode,thisNode,thisWeight);
  180.   thisNode->parentEdge = thisEdge;
  181.   if (NULL == centerNode->leftEdge)
  182.     centerNode->leftEdge = thisEdge;
  183.   else if (NULL == centerNode->rightEdge)
  184.     centerNode->rightEdge = thisEdge;
  185.   else if (NULL == centerNode->middleEdge)
  186.     centerNode->middleEdge = thisEdge;
  187.   else
  188.     {
  189.       fprintf(stderr,"Error: meNode %s has too many (>3) children.n",centerNode->label);
  190.       exit(EXIT_FAILURE);
  191.     }
  192.   sprintf(thisEdge->label,"E%d",edgeCount++);
  193.   i = right + 1;
  194.   break;
  195. }
  196.     }
  197.   return(centerNode);
  198. }
  199. meTree *loadNewickTree(FILE *ifile, int numLeaves)
  200. {
  201.   char label[] = "EmptyEdge";
  202.   meTree *T;
  203.   meNode *centerNode;
  204.   int i = 0;
  205.   int j = 0;
  206.   int inputLength;
  207.   int uCount = 0;
  208.   int parCount = 0;
  209.   char c;
  210.   boolean Comment;
  211.   char *nextString;
  212.   char rootLabel[MAX_LABEL_LENGTH];
  213.   nodeCount = edgeCount = 0;
  214.   T = newTree();
  215.   nextString = (char *) malloc(numLeaves*INPUT_SIZE*sizeof(char));
  216.   if (NULL == nextString)
  217.     nextString = (char *) malloc(MAX_INPUT_SIZE*sizeof(char));
  218.   Comment = FALSE_FASTME;
  219.   while(1 == fscanf(ifile,"%c",&c))
  220.     {
  221.       if('[' == c)
  222. Comment = TRUE_FASTME;
  223.       else if (']' == c)
  224. Comment = FALSE_FASTME;
  225.       else if (!(Comment))
  226. {
  227.   if(whiteSpace(c)) 
  228.     {
  229.       if (i > 0)
  230. nextString[i++] = ' ';
  231.     }
  232.   else  /*note that this else goes with if(whiteSpace(c))*/
  233.     nextString[i++] = c;
  234.   if (';' == c)
  235.     break;
  236. }
  237.     }
  238.   if ('(' != nextString[0])
  239.     {
  240.       fprintf(stderr,"Error reading input file - does not start with '('.n");
  241.       exit(EXIT_FAILURE);
  242.     }
  243.   inputLength = i;
  244.   for(i = 0; i < inputLength;i++)
  245.     {
  246.       if ('(' == nextString[i])
  247. parCount++;
  248.       else if (')' == nextString[i])
  249. parCount--;
  250.       if (parCount > 0)
  251. ;
  252.       else if (0 == parCount)
  253. {
  254.   i++;
  255.   if(';' == nextString[i])
  256.     sprintf(rootLabel,"URoot");  
  257.   else
  258.     {
  259.       while(';' != nextString[i]) 
  260. if(!(whiteSpace(nextString[i++])))
  261.   rootLabel[j++] = nextString[i-1];  /*be careful here */
  262.       rootLabel[j] = '';
  263.     }
  264.   i = inputLength;
  265. }
  266.       else if (parCount < 0)
  267. {
  268.   fprintf(stderr,"Error reading meTree input file.  Too many right parentheses.n");
  269.   exit(EXIT_FAILURE);
  270. }
  271.     }
  272.   centerNode = decodeNewickSubtree(nextString,T,&uCount);
  273.   sprintf(centerNode->label,rootLabel);
  274.   T->root = centerNode;
  275.   free(nextString);
  276.   return(T);
  277. }
  278. void NewickPrintSubtree(meTree *T, meEdge *e, FILE *ofile)
  279. {
  280.   if (NULL == e)
  281.     {
  282.       fprintf(stderr,"Error with Newick Printing routine.n");
  283.       exit(EXIT_FAILURE);
  284.     }
  285.   if(!(leaf(e->head)))
  286.     {
  287.       fprintf(ofile,"(");
  288.       NewickPrintSubtree(T,e->head->leftEdge,ofile);
  289.       fprintf(ofile,", ");
  290.       NewickPrintSubtree(T,e->head->rightEdge,ofile);
  291.       fprintf(ofile,")");
  292.     }
  293.   fprintf(ofile," ");
  294.   fprintf(ofile,"%s ",e->head->label);
  295.   fprintf(ofile,":%lf",e->distance);
  296. }
  297.       
  298. int NewickPrintSubtreeString(meTree *T, meEdge *e, char *s)
  299. {
  300.   meEdge *left, *right;
  301.   int slength, llength, rlength;
  302.   slength = 0;
  303.   left = e->head->leftEdge;
  304.   right = e->head->rightEdge;
  305.   if (!(NULL == left))
  306.     {
  307.       s[0] = '(';
  308.       left = e->head->leftEdge;
  309.       llength = NewickPrintSubtreeString(T,left,s+1);
  310.       s[1 + llength] = ',';
  311.       right = e->head->rightEdge;
  312.       rlength = NewickPrintSubtreeString(T,right,s+2+llength);
  313.       s[2 + rlength + llength]=')';
  314.       slength = 3 + rlength + llength;
  315.     }
  316.   sprintf(s+slength,"%s:%lf",e->head->label,e->distance);
  317.   while ('' != s[slength])
  318.     slength++;
  319.   return(slength);
  320. }
  321. void NewickPrintTree(meTree *T, FILE *ofile)
  322. {
  323.   meEdge *e, *f;
  324.   meNode *rootchild;
  325.   e = T->root->leftEdge;
  326.   rootchild = e->head;
  327.   fprintf(ofile,"(%s: %lf",T->root->label,e->distance);
  328.   f = rootchild->leftEdge;
  329.   if (NULL != f)
  330.     {
  331.       fprintf(ofile,",");
  332.       NewickPrintSubtree(T,f,ofile);
  333.     }
  334.   f = rootchild->rightEdge;
  335.   if (NULL != f)
  336.     {
  337.       fprintf(ofile,",");
  338.       NewickPrintSubtree(T,f,ofile);
  339.     }
  340.   fprintf(ofile,")");
  341.   if (NULL != rootchild->label)
  342.     fprintf(ofile,"%s",rootchild->label);
  343.   fprintf(ofile,";n");
  344. }
  345. char *NewickPrintTreeString(meTree *T)
  346. {
  347.   char *s;
  348.   int i,slength, llength, rlength;
  349.   meEdge *e, *left, *right;
  350.   e = T->root->leftEdge;
  351.   slength = T->size * NODE_LABEL_LENGTH + (T->size-1)*MAX_DIGITS
  352.     + 2*T->size + T->size + (T->size-1) + 1;
  353.   /*sum represents characters needed for T->size meNode labels, T->size-1 
  354.     meEdge lengths, at most 2*T->size parentheses, at most T->size commas, 
  355.     T->size - 1 ':' characters for meEdge lengths, and 1 ';' character at end*/
  356.   
  357.   s = (char *) malloc(slength*sizeof(char));
  358.   left = e->head->leftEdge;
  359.   right = e->head->rightEdge;
  360.   i=0;
  361.   s[0] = '(';
  362.   if (!(NULL == left))
  363.     {
  364.       left = e->head->leftEdge;
  365.       llength = NewickPrintSubtreeString(T,left,s+1);
  366.       s[1 + llength] = ',';
  367.       right = e->head->rightEdge;
  368.       rlength = NewickPrintSubtreeString(T,right,s+2+llength);
  369.       s[2 + rlength + llength]=',';
  370.       i = 3 + rlength + llength;
  371.     }
  372.   sprintf(s+i,"%s:%lf)%s;",e->tail->label,e->distance,e->head->label);
  373.   return(s);
  374. }
  375. meEdge *depthFirstTraverse(meTree *T, meEdge *e);
  376. END_SCOPE(fastme)
  377. END_NCBI_SCOPE
  378. /*
  379.  * ===========================================================================
  380.  * $Log: newickstring.cpp,v $
  381.  * Revision 1000.1  2004/06/01 18:10:14  gouriano
  382.  * PRODUCTION: UPGRADED [GCC34_MSVC7] Dev-tree R1.2
  383.  *
  384.  * Revision 1.2  2004/05/21 21:41:04  gorelenk
  385.  * Added PCH ncbi_pch.hpp
  386.  *
  387.  * Revision 1.1  2004/02/10 15:16:03  jcherry
  388.  * Initial version
  389.  *
  390.  * ===========================================================================
  391.  */