5 * Created by Pat Schloss on 6/17/09.
6 * Copyright 2009 Patrick D. Schloss. All rights reserved.
10 #include "doTaxonomy.h"
12 /**************************************************************************************************/
14 PhyloTree::PhyloTree(){
18 tree.push_back(TaxNode("Root"));
19 tree[0].heirarchyID = "0";
22 /**************************************************************************************************/
24 string PhyloTree::getNextTaxon(string& heirarchy){
26 string currentLevel = "";
28 currentLevel=heirarchy.substr(0,heirarchy.find_first_of(';'));
29 heirarchy=heirarchy.substr(heirarchy.find_first_of(';')+1);
35 /**************************************************************************************************/
37 void PhyloTree::addSeqToTree(string seqName, string seqTaxonomy){
41 map<string, int>::iterator childPointer;
46 tree[0].accessions.push_back(seqName);
47 string taxon;// = getNextTaxon(seqTaxonomy);
49 while(seqTaxonomy != ""){
53 //somehow the parent is getting one too many accnos
54 //use print to reassign the taxa id
55 taxon = getNextTaxon(seqTaxonomy);
57 childPointer = tree[currentNode].children.find(taxon);
59 if(childPointer != tree[currentNode].children.end()){ //if the node already exists, move on
60 currentNode = childPointer->second;
61 tree[currentNode].accessions.push_back(seqName);
62 name2Taxonomy[seqName] = currentNode;
64 else{ //otherwise, create it
65 tree.push_back(TaxNode(taxon));
67 tree[currentNode].children[taxon] = numNodes-1;
68 tree[numNodes-1].parent = currentNode;
70 // int numChildren = tree[currentNode].children.size();
71 // string heirarchyID = tree[currentNode].heirarchyID;
72 // tree[currentNode].accessions.push_back(seqName);
74 currentNode = tree[currentNode].children[taxon];
75 tree[currentNode].accessions.push_back(seqName);
76 name2Taxonomy[seqName] = currentNode;
77 // tree[currentNode].level = level;
78 // tree[currentNode].childNumber = numChildren;
79 // tree[currentNode].heirarchyID = heirarchyID + '.' + toString(tree[currentNode].childNumber);
82 if (seqTaxonomy == "") { uniqueTaxonomies[currentNode] = currentNode; }
85 /**************************************************************************************************/
86 vector<int> PhyloTree::getGenusNodes() {
89 //generate genusIndexes
90 map<int, int>::iterator it2;
91 for (it2=uniqueTaxonomies.begin(); it2!=uniqueTaxonomies.end(); it2++) { genusIndex.push_back(it2->first); }
96 errorOut(e, "PhyloTree", "getGenusNodes");
100 /**************************************************************************************************/
102 void PhyloTree::assignHeirarchyIDs(int index){
104 map<string,int>::iterator it;
107 for(it=tree[index].children.begin();it!=tree[index].children.end();it++){
108 tree[it->second].heirarchyID = tree[index].heirarchyID + '.' + toString(counter);
110 tree[it->second].level = tree[index].level + 1;
111 assignHeirarchyIDs(it->second);
115 /**************************************************************************************************/
117 void PhyloTree::print(ofstream& out){
119 out << tree[0].level << '\t'<< tree[0].heirarchyID << '\t' << tree[0].name << '\t' << tree[0].children.size() << '\t' << tree[0].accessions.size() << endl;
125 /**************************************************************************************************/
127 void PhyloTree::print(int i, ofstream& out){
129 map<string,int>::iterator it;
130 for(it=tree[i].children.begin();it!=tree[i].children.end();it++){
131 out <<tree[it->second].level << '\t' << tree[it->second].heirarchyID << '\t' << tree[it->second].name << '\t' << tree[it->second].children.size() << '\t' << tree[it->second].accessions.size() << endl;
132 print(it->second, out);
137 /**************************************************************************************************/