From: westcott Date: Thu, 5 Mar 2009 19:54:36 +0000 (+0000) Subject: working on readtree X-Git-Url: https://git.donarmstrong.com/?p=mothur.git;a=commitdiff_plain;h=8f92612d5a69f5245e63a20657e7d93519c0769e working on readtree --- diff --git a/clustercommand.cpp b/clustercommand.cpp index dc61a6a..f0ef51f 100644 --- a/clustercommand.cpp +++ b/clustercommand.cpp @@ -106,7 +106,14 @@ int ClusterCommand::execute(){ else if(rndPreviousDistsetSparseMatrix(clearM); + ListVector* clearL = NULL; + globaldata->setListVector(clearL); + + //saves .list file so you can do the collect, rarefaction and summary commands without doing a read.list if (globaldata->getFormat() == "phylip") { globaldata->setPhylipFile(""); } else if (globaldata->getFormat() == "column") { globaldata->setColumnFile(""); } diff --git a/errorchecking.cpp b/errorchecking.cpp index ff94c8f..b76542a 100644 --- a/errorchecking.cpp +++ b/errorchecking.cpp @@ -367,8 +367,8 @@ bool ErrorCheck::checkInput(string input) { } //are you trying to cluster before you have read something - if ((commandName == "cluster") && (globaldata->getSparseMatrix() == NULL) || - (commandName == "cluster") && (globaldata->getListVector() == NULL)) { + if (((commandName == "cluster") && (globaldata->getSparseMatrix() == NULL)) || + ((commandName == "cluster") && (globaldata->getListVector() == NULL))) { cout << "Before you use the cluster command, you first need to read in a distance matrix." << endl; errorFree = false; } diff --git a/readtree.cpp b/readtree.cpp index 81a7729..4d3bd6d 100644 --- a/readtree.cpp +++ b/readtree.cpp @@ -124,7 +124,8 @@ void ReadNewickTree::read() { numNodes = T->getNumNodes(); numLeaves = T->getNumLeaves(); - readTreeString(); + readTreeString(); + //save trees for later commands globaldata->gTree.push_back(T); gobble(filehandle); @@ -157,12 +158,12 @@ void ReadNewickTree::read() { numLeaves = T->getNumLeaves(); //read tree info - readTreeString(); + readTreeString(); + //save trees for later commands globaldata->gTree.push_back(T); } } - } catch(exception& e) { cout << "Standard Error: " << e.what() << " has occurred in the ReadNewickTree class Function read. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n"; @@ -272,6 +273,7 @@ void ReadNewickTree::readTreeString() { T->tree[n].setParent(-1); if(lc!=-1){ T->tree[lc].setParent(n); } if(rc!=-1){ T->tree[rc].setParent(n); } + cout << "new loop "<< endl; } } @@ -334,6 +336,7 @@ int ReadNewickTree::readNewickInt(istream& f, int& n, Tree* T) { //adds sequence names that are not in group file to the "xxx" group if(n1 == -1) { cerr << "Name: " << name << " not found in your groupfile and it will be ignored. \n"; + globaldata->gTreemap->namesOfSeqs.push_back(name); globaldata->gTreemap->treemap[name].groupname = "xxx"; globaldata->gTreemap->treemap[name].vectorIndex = (globaldata->gTreemap->namesOfSeqs.size() - 1); @@ -350,10 +353,13 @@ int ReadNewickTree::readNewickInt(istream& f, int& n, Tree* T) { //find index in tree of name n1 = T->getIndex(name); group = "xxx"; + numLeaves++; + numNodes = 2*numLeaves - 1; + T->resetTree(); } T->tree[n1].setGroup(group); - + T->printTree(); T->tree[n1].setChildren(-1,-1); if(blen == 1){ diff --git a/readtree.h b/readtree.h index bee4d94..250b1db 100644 --- a/readtree.h +++ b/readtree.h @@ -43,7 +43,7 @@ class ReadTree { class ReadNewickTree : public ReadTree { public: - ReadNewickTree(string file) : treeFile(file) { openInputFile(file, filehandle); } + ReadNewickTree(string file) : treeFile(file) { openInputFile(file, filehandle); } ~ReadNewickTree() {}; void read(); diff --git a/readtreecommand.cpp b/readtreecommand.cpp index a540c74..68f1f92 100644 --- a/readtreecommand.cpp +++ b/readtreecommand.cpp @@ -52,7 +52,7 @@ int ReadTreeCommand::execute(){ for (int i = 0; i < T.size(); i++) { T[i]->assembleTree(); } - +T[0]->createNewickFile("treeout"); return 0; } catch(exception& e) { diff --git a/tree.cpp b/tree.cpp index 642c665..e34253d 100644 --- a/tree.cpp +++ b/tree.cpp @@ -26,8 +26,6 @@ Tree::Tree() { if (i <= (numLeaves-1)) { tree[i].setName(globaldata->gTreemap->namesOfSeqs[i]); tree[i].setGroup(globaldata->gTreemap->getGroup(globaldata->gTreemap->namesOfSeqs[i])); - //the node knows its index - tree[i].setIndex(i); //set pcount and pGroup for groupname to 1. tree[i].pcount[globaldata->gTreemap->getGroup(globaldata->gTreemap->namesOfSeqs[i])] = 1; tree[i].pGroups[globaldata->gTreemap->getGroup(globaldata->gTreemap->namesOfSeqs[i])] = 1; @@ -38,8 +36,6 @@ Tree::Tree() { }else if (i > (numLeaves-1)) { tree[i].setName(""); tree[i].setGroup(""); - //the node knows its index - tree[i].setIndex(i); } } } @@ -54,7 +50,43 @@ Tree::Tree() { } /*****************************************************************/ +void Tree::resetTree(){ +try { + numLeaves = globaldata->gTreemap->getNumSeqs(); + numNodes = 2*numLeaves - 1; + + tree.resize(numNodes); + //initialize tree with correct number of nodes, name and group info. + for (int i = 0; i < numNodes; i++) { + //initialize leaf nodes + if (i <= (numLeaves-1)) { + tree[i].setName(globaldata->gTreemap->namesOfSeqs[i]); + tree[i].setGroup(globaldata->gTreemap->getGroup(globaldata->gTreemap->namesOfSeqs[i])); + //set pcount and pGroup for groupname to 1. + tree[i].pcount[globaldata->gTreemap->getGroup(globaldata->gTreemap->namesOfSeqs[i])] = 1; + tree[i].pGroups[globaldata->gTreemap->getGroup(globaldata->gTreemap->namesOfSeqs[i])] = 1; + //Treemap knows name, group and index to speed up search + globaldata->gTreemap->setIndex(globaldata->gTreemap->namesOfSeqs[i], i); + + //intialize non leaf nodes + }else if (i > (numLeaves-1)) { + tree[i].setName(""); + tree[i].setGroup(""); + } + } + } + catch(exception& e) { + cout << "Standard Error: " << e.what() << " has occurred in the Tree class Function resetTree. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n"; + exit(1); + } + catch(...) { + cout << "An unknown error has occurred in the Tree class function resetTree. Please contact Pat Schloss at pschloss@microbio.umass.edu." << "\n"; + exit(1); + } + +} +/*****************************************************************/ int Tree::getIndex(string searchName) { try { //Treemap knows name, group and index to speed up search @@ -157,7 +189,7 @@ map Tree::mergeGroups(int i) { try { int lc = tree[i].getLChild(); int rc = tree[i].getRChild(); - +cout << i << lc << rc << endl; //set parsimony groups to left child map parsimony = tree[lc].pGroups; diff --git a/tree.h b/tree.h index d6a3b55..828b038 100644 --- a/tree.h +++ b/tree.h @@ -24,6 +24,7 @@ public: void getCopy(Tree*); //makes tree a copy of the one passed in. + void resetTree(); //this is needed to allow user to omit names from the group file void assembleRandomTree(); void assembleRandomUnifracTree(vector); void assembleRandomUnifracTree(string, string); diff --git a/treemap.cpp b/treemap.cpp index 6314eb0..9eb8330 100644 --- a/treemap.cpp +++ b/treemap.cpp @@ -47,7 +47,7 @@ void TreeMap::readMap() { int TreeMap::getNumGroups() { - return seqsPerGroup.size(); + return namesOfGroups.size(); } /************************************************************/ diff --git a/unifracunweightedcommand.cpp b/unifracunweightedcommand.cpp index 2a169d0..3fa44b5 100644 --- a/unifracunweightedcommand.cpp +++ b/unifracunweightedcommand.cpp @@ -22,7 +22,6 @@ UnifracUnweightedCommand::UnifracUnweightedCommand() { setGroups(); //sets users groups to analyze convert(globaldata->getIters(), iters); //how many random trees to generate unweighted = new Unweighted(tmap); - counter = 0; } catch(exception& e) { @@ -44,6 +43,7 @@ int UnifracUnweightedCommand::execute() { //get pscores for users trees for (int i = 0; i < T.size(); i++) { + counter = 0; unweightedFile = globaldata->getTreeFile() + toString(i+1) + ".unweighted"; unweightedFileout = globaldata->getTreeFile() + "temp." + toString(i+1) + ".unweighted"; //column headers @@ -53,7 +53,6 @@ int UnifracUnweightedCommand::execute() { //get unweighted for users tree rscoreFreq.resize(numComp); rCumul.resize(numComp); - validScores.resize(numComp); utreeScores.resize(numComp); UWScoreSig.resize(numComp); @@ -62,9 +61,6 @@ int UnifracUnweightedCommand::execute() { //output scores for each combination for(int k = 0; k < numComp; k++) { - //add users score to valid scores - validScores[k][userData[k]] = userData[k]; - //saves users score utreeScores[k].push_back(userData[k]); } @@ -85,14 +81,14 @@ int UnifracUnweightedCommand::execute() { } //add randoms score to validscores - validScores[k][randomData[k]] = randomData[k]; + validScores[randomData[k]] = randomData[k]; } } for(int a = 0; a < numComp; a++) { float rcumul = 1.0000; //this loop fills the cumulative maps and put 0.0000 in the score freq map to make it easier to print. - for (it = validScores[a].begin(); it != validScores[a].end(); it++) { + for (it = validScores.begin(); it != validScores.end(); it++) { //make rscoreFreq map and rCumul it2 = rscoreFreq[a].find(it->first); rCumul[a][it->first] = rcumul; @@ -136,7 +132,7 @@ void UnifracUnweightedCommand::printUnweightedFile() { for(int a = 0; a < numComp; a++) { initFile(groupComb[a]); //print each line - for (it = validScores[a].begin(); it != validScores[a].end(); it++) { + for (it = validScores.begin(); it != validScores.end(); it++) { data.push_back(it->first); data.push_back(rscoreFreq[a][it->first]); data.push_back(rCumul[a][it->first]); output(data); data.clear();