]> git.donarmstrong.com Git - mothur.git/blobdiff - splitmatrix.cpp
working on adding count file to mgcluster
[mothur.git] / splitmatrix.cpp
index cd05cc46c172122ab9ee23490d38eda5e7308c22..384b09af1bb94be09c5607c8b863b28b16215731 100644 (file)
@@ -9,9 +9,8 @@
 
 #include "splitmatrix.h"
 #include "phylotree.h"
-#include "sequencedb.h"
-#include "onegapdist.h"
-#include "dist.h"
+#include "distancecommand.h"
+#include "seqsummarycommand.h"
 
 /***********************************************************************/
 
@@ -26,12 +25,17 @@ SplitMatrix::SplitMatrix(string distfile, string name, string tax, float c, stri
 }
 /***********************************************************************/
 
-SplitMatrix::SplitMatrix(string ffile, string tax, float c, string t){
+SplitMatrix::SplitMatrix(string ffile, string name, string tax, float c, float cu, string t, int p, bool cl, string output){
        m = MothurOut::getInstance();
        fastafile = ffile;
+       namefile = name;
        taxFile = tax;
-       cutoff = c;
+       cutoff = c;  //tax level cutoff
+       distCutoff = cu; //for fasta method if you are creating distance matrix you need a cutoff for that
        method = t;
+       processors = p;
+    classic = cl;
+       outputDir = output;
 }
 
 /***********************************************************************/
@@ -63,6 +67,8 @@ int SplitMatrix::splitDistance(){
         
                if (large)      { splitDistanceLarge(); }
                else            { splitDistanceRAM();   }
+               
+               return 0;
                        
        }
        catch(exception& e) {
@@ -75,7 +81,7 @@ int SplitMatrix::splitDistance(){
 int SplitMatrix::splitClassify(){
        try {
                cutoff = int(cutoff);
-               
+                               
                map<string, int> seqGroup;
                map<string, int>::iterator it;
                map<string, int>::iterator it2;
@@ -85,16 +91,13 @@ int SplitMatrix::splitClassify(){
                //build tree from users taxonomy file
                PhyloTree* phylo = new PhyloTree();
                
-               ifstream in;
-               openInputFile(taxFile, in);
-                       
-               //read in users taxonomy file and add sequences to tree
-               string seqname, tax;
-               while(!in.eof()){
-                       in >> seqname >> tax; gobble(in);
-                       phylo->addSeqToTree(seqname, tax);
-               }
-               in.close();
+        map<string, string> temp;
+        m->readTax(taxFile, temp);
+        
+        for (map<string, string>::iterator itTemp = temp.begin(); itTemp != temp.end();) {
+            phylo->addSeqToTree(itTemp->first, itTemp->second);
+            temp.erase(itTemp++);
+        }
                
                phylo->assignHeirarchyIDs(0);
 
@@ -125,7 +128,6 @@ int SplitMatrix::splitClassify(){
                        createDistanceFilesFromTax(seqGroup, numGroups);
                }
                
-                               
                return 0;
                        
        }
@@ -137,36 +139,132 @@ int SplitMatrix::splitClassify(){
 /***********************************************************************/
 int SplitMatrix::createDistanceFilesFromTax(map<string, int>& seqGroup, int numGroups){
        try {
+               map<string, int> copyGroups = seqGroup;
                map<string, int>::iterator it;
-               map<string, int>::iterator it2;
-               map<string, int> seqIndexInFasta;
+               set<string> names;
                                
-               //read fastafile
-               SequenceDB alignDB;
-               
-               ifstream filehandle;
-               openInputFile(fastafile, filehandle);
-               int numSeqs = 0;
-               while (!filehandle.eof()) {
-                       //input sequence info into sequencedb
-                       Sequence newSequence(filehandle);
+               for (int i = 0; i < numGroups; i++) { //remove old temp files, just in case
+                       m->mothurRemove((fastafile + "." + toString(i) + ".temp"));
+               }
                        
-                       if (newSequence.getName() != "") {   
-                               alignDB.push_back(newSequence);  
-                               seqIndexInFasta[newSequence.getName()] = numSeqs;
-                               numSeqs++;
+               ifstream in;
+               m->openInputFile(fastafile, in);
+       
+               //parse fastafile
+               ofstream outFile;
+               while (!in.eof()) {
+                       Sequence query(in); m->gobble(in);
+                       if (query.getName() != "") {
+               
+                               it = seqGroup.find(query.getName());
+                               
+                               //save names in case no namefile is given
+                               if (namefile == "") {  names.insert(query.getName()); }
+                       
+                               if (it != seqGroup.end()) { //not singleton 
+                                       m->openOutputFileAppend((fastafile + "." + toString(it->second) + ".temp"), outFile);
+                                       query.printSequence(outFile); 
+                                       outFile.close();
+                                       
+                                       copyGroups.erase(query.getName());
+                               }
                        }
+               }
+               in.close();
+               
+               //warn about sequence in groups that are not in fasta file
+               for(it = copyGroups.begin(); it != copyGroups.end(); it++) {
+                       m->mothurOut("ERROR: " + it->first + " is missing from your fastafile. This could happen if your taxonomy file is not unique and your fastafile is, or it could indicate and error."); m->mothurOutEndLine();
+                       exit(1);
+               }
+               
+               copyGroups.clear();
+               
+               //process each distance file
+               for (int i = 0; i < numGroups; i++) { 
+                       
+                       string options = "";
+            if (classic) { options = "fasta=" + (fastafile + "." + toString(i) + ".temp") + ", processors=" + toString(processors) + ", output=lt"; }
+            else { options = "fasta=" + (fastafile + "." + toString(i) + ".temp") + ", processors=" + toString(processors) + ", cutoff=" + toString(distCutoff); }
+                       if (outputDir != "") { options += ", outputdir=" + outputDir; }
+                       
+                       Command* command = new DistanceCommand(options);
+                       
+                       command->execute();
+                       delete command;
                        
-                       //takes care of white space
-                       gobble(filehandle);
+                       m->mothurRemove((fastafile + "." + toString(i) + ".temp"));
+                       
+                       //remove old names files just in case
+                       m->mothurRemove((namefile + "." + toString(i) + ".temp"));
                }
-               filehandle.close();
+                       
+               singleton = namefile + ".extra.temp";
+               ofstream remainingNames;
+               m->openOutputFile(singleton, remainingNames);
                
-               Dist* distCalculator = new oneGapDist();
+               bool wroteExtra = false;
+
+               ifstream bigNameFile;
+               m->openInputFile(namefile, bigNameFile);
                
+               string name, nameList;
+               while(!bigNameFile.eof()){
+                       bigNameFile >> name >> nameList;  m->gobble(bigNameFile);
+                       
+                       //did this sequence get assigned a group
+                       it = seqGroup.find(name);
+                       
+                       if (it != seqGroup.end()) {  
+                               m->openOutputFileAppend((namefile + "." + toString(it->second) + ".temp"), outFile);
+                               outFile << name << '\t' << nameList << endl;
+                               outFile.close();
+                       }else{
+                               wroteExtra = true;
+                               remainingNames << name << '\t' << nameList << endl;
+                       }
+               }
+               bigNameFile.close();
                
-//still not done....           
+               for(int i=0;i<numGroups;i++){
+                       string tempNameFile = namefile + "." + toString(i) + ".temp";
+                       if (outputDir == "") { outputDir = m->hasPath(fastafile); }
+                       string tempDistFile = "";
+            if (classic) { tempDistFile =  outputDir + m->getRootName(m->getSimpleName((fastafile + "." + toString(i) + ".temp"))) + "phylip.dist";}
+            else { tempDistFile = outputDir + m->getRootName(m->getSimpleName((fastafile + "." + toString(i) + ".temp"))) + "dist"; }
+
+                       //if there are valid distances
+                       ifstream fileHandle;
+                       fileHandle.open(tempDistFile.c_str());
+                       if(fileHandle)  {       
+                               m->gobble(fileHandle);
+                               if (!fileHandle.eof()) {  //check for blank file - this could occur if all dists in group are above cutoff
+                                       map<string, string> temp;
+                                       temp[tempDistFile] = tempNameFile;
+                                       dists.push_back(temp);
+                               }else {
+                                       ifstream in;
+                                       m->openInputFile(tempNameFile, in);
+                               
+                                       while(!in.eof()) { 
+                                               in >> name >> nameList;  m->gobble(in);
+                                               wroteExtra = true;
+                                               remainingNames << name << '\t' << nameList << endl;
+                                       }
+                                       in.close();
+                                       m->mothurRemove(tempNameFile);
+                               }
+                       }
+                       fileHandle.close();
+               }
                
+               remainingNames.close();
+               if (!wroteExtra) { 
+                       m->mothurRemove(singleton);
+                       singleton = "none";
+               }
+
+               if (m->control_pressed)  {  for (int i = 0; i < dists.size(); i++) { m->mothurRemove((dists[i].begin()->first)); m->mothurRemove((dists[i].begin()->second)); } dists.clear(); }
                
                return 0;
        }
@@ -182,11 +280,11 @@ int SplitMatrix::splitDistanceFileByTax(map<string, int>& seqGroup, int numGroup
                map<string, int>::iterator it2;
                
                ifstream dFile;
-               openInputFile(distFile, dFile);
+               m->openInputFile(distFile, dFile);
                ofstream outFile;
                
                for (int i = 0; i < numGroups; i++) { //remove old temp files, just in case
-                       remove((distFile + "." + toString(i) + ".temp").c_str());
+                       m->mothurRemove((distFile + "." + toString(i) + ".temp"));
                }
                
                //for buffering the io to improve speed
@@ -203,9 +301,9 @@ int SplitMatrix::splitDistanceFileByTax(map<string, int>& seqGroup, int numGroup
                        string seqA, seqB;
                        float dist;
                        
-                       if (m->control_pressed) { dFile.close(); for (int i = 0; i < numGroups; i++) { remove((distFile + "." + toString(i) + ".temp").c_str());        } }
+                       if (m->control_pressed) { dFile.close(); for (int i = 0; i < numGroups; i++) { m->mothurRemove((distFile + "." + toString(i) + ".temp"));       } }
                        
-                       dFile >> seqA >> seqB >> dist;  gobble(dFile);
+                       dFile >> seqA >> seqB >> dist;  m->gobble(dFile);
                        
                        //if both sequences are in the same group then they are within the cutoff
                        it = seqGroup.find(seqA);
@@ -214,7 +312,7 @@ int SplitMatrix::splitDistanceFileByTax(map<string, int>& seqGroup, int numGroup
                        if ((it != seqGroup.end()) && (it2 != seqGroup.end())) { //they are both not singletons 
                                if (it->second == it2->second) { //they are from the same group so add the distance
                                        if (numOutputs[it->second] > 30) {
-                                               openOutputFileAppend((distFile + "." + toString(it->second) + ".temp"), outFile);
+                                               m->openOutputFileAppend((distFile + "." + toString(it->second) + ".temp"), outFile);
                                                outFile << outputs[it->second] << seqA << '\t' << seqB << '\t' << dist << endl;
                                                outFile.close();
                                                outputs[it->second] = "";
@@ -230,11 +328,11 @@ int SplitMatrix::splitDistanceFileByTax(map<string, int>& seqGroup, int numGroup
                dFile.close();
        
                for (int i = 0; i < numGroups; i++) { //remove old temp files, just in case
-                       remove((namefile + "." + toString(i) + ".temp").c_str());
+                       m->mothurRemove((namefile + "." + toString(i) + ".temp"));
                        
                        //write out any remaining buffers
                        if (numOutputs[i] > 0) {
-                               openOutputFileAppend((distFile + "." + toString(i) + ".temp"), outFile);
+                               m->openOutputFileAppend((distFile + "." + toString(i) + ".temp"), outFile);
                                outFile << outputs[i];
                                outFile.close();
                                outputs[i] = "";
@@ -244,23 +342,23 @@ int SplitMatrix::splitDistanceFileByTax(map<string, int>& seqGroup, int numGroup
                }
                
                ifstream bigNameFile;
-               openInputFile(namefile, bigNameFile);
+               m->openInputFile(namefile, bigNameFile);
                
                singleton = namefile + ".extra.temp";
                ofstream remainingNames;
-               openOutputFile(singleton, remainingNames);
+               m->openOutputFile(singleton, remainingNames);
                
                bool wroteExtra = false;
                                                
                string name, nameList;
                while(!bigNameFile.eof()){
-                       bigNameFile >> name >> nameList;  gobble(bigNameFile);
+                       bigNameFile >> name >> nameList;  m->gobble(bigNameFile);
                        
                        //did this sequence get assigned a group
                        it = seqGroup.find(name);
                        
                        if (it != seqGroup.end()) {  
-                               openOutputFileAppend((namefile + "." + toString(it->second) + ".temp"), outFile);
+                               m->openOutputFileAppend((namefile + "." + toString(it->second) + ".temp"), outFile);
                                outFile << name << '\t' << nameList << endl;
                                outFile.close();
                        }else{
@@ -269,29 +367,41 @@ int SplitMatrix::splitDistanceFileByTax(map<string, int>& seqGroup, int numGroup
                        }
                }
                bigNameFile.close();
-               remainingNames.close();
-               
-               if (!wroteExtra) { 
-                       remove(singleton.c_str());
-                       singleton = "none";
-               }
-               
+                               
                for(int i=0;i<numGroups;i++){
+                       string tempNameFile = namefile + "." + toString(i) + ".temp";
+                       string tempDistFile = distFile + "." + toString(i) + ".temp";
+
                        //if there are valid distances
                        if (validDistances[i]) {
-                               string tempNameFile = namefile + "." + toString(i) + ".temp";
-                               string tempDistFile = distFile + "." + toString(i) + ".temp";
-                               
                                map<string, string> temp;
                                temp[tempDistFile] = tempNameFile;
                                dists.push_back(temp);
+                       }else{
+                               ifstream in;
+                               m->openInputFile(tempNameFile, in);
+                               
+                               while(!in.eof()) { 
+                                       in >> name >> nameList;  m->gobble(in);
+                                       wroteExtra = true;
+                                       remainingNames << name << '\t' << nameList << endl;
+                               }
+                               in.close();
+                               m->mothurRemove(tempNameFile);
                        }
                }
                
+               remainingNames.close();
+               
+               if (!wroteExtra) { 
+                       m->mothurRemove(singleton);
+                       singleton = "none";
+               }
+
                if (m->control_pressed)  {  
                        for (int i = 0; i < dists.size(); i++) { 
-                               remove((dists[i].begin()->first).c_str());
-                               remove((dists[i].begin()->second).c_str());
+                               m->mothurRemove((dists[i].begin()->first));
+                               m->mothurRemove((dists[i].begin()->second));
                        }
                        dists.clear();
                }
@@ -318,7 +428,7 @@ int SplitMatrix::splitDistanceLarge(){
 
                ofstream outFile;
                ifstream dFile;
-               openInputFile(distFile, dFile);
+               m->openInputFile(distFile, dFile);
        
                while(dFile){
                        string seqA, seqB;
@@ -326,7 +436,7 @@ int SplitMatrix::splitDistanceLarge(){
 
                        dFile >> seqA >> seqB >> dist;
                        
-                       if (m->control_pressed) {   dFile.close();  for(int i=0;i<numGroups;i++){       if(groups[i].size() > 0){  remove((distFile + "." + toString(i) + ".temp").c_str()); }  } return 0; }
+                       if (m->control_pressed) {   dFile.close();  for(int i=0;i<numGroups;i++){       if(groups[i].size() > 0){  m->mothurRemove((distFile + "." + toString(i) + ".temp")); }  } return 0; }
                                        
                        if(dist < cutoff){
                                //cout << "in cutoff: " << dist << endl;
@@ -454,7 +564,7 @@ int SplitMatrix::splitDistanceLarge(){
                                                                delete memblock;
                                                                
                                                                fileB.close();
-                                                               remove(fileName2.c_str());
+                                                               m->mothurRemove(fileName2);
                                                                
                                                                //write out the merged memory
                                                                if (numOutputs[groupID] > 60) {
@@ -514,7 +624,7 @@ int SplitMatrix::splitDistanceLarge(){
                                                                delete memblock;
                                                                
                                                                fileB.close();
-                                                               remove(fileName2.c_str());
+                                                               m->mothurRemove(fileName2);
                                                                
                                                                //write out the merged memory
                                                                if (numOutputs[groupID] > 60) {
@@ -532,7 +642,7 @@ int SplitMatrix::splitDistanceLarge(){
                                        }
                                }
                        }
-                       gobble(dFile);
+                       m->gobble(dFile);
                }
                dFile.close();
                
@@ -570,7 +680,7 @@ int SplitMatrix::splitNames(vector<set<string> >& groups){
                while(bigNameFile){
                        bigNameFile >> name >> nameList;
                        nameMap[name] = nameList;
-                       gobble(bigNameFile);
+                       m->gobble(bigNameFile);
                }
                bigNameFile.close();
                        
@@ -617,8 +727,8 @@ int SplitMatrix::splitNames(vector<set<string> >& groups){
                
                if (m->control_pressed)  {  
                        for (int i = 0; i < dists.size(); i++) { 
-                               remove((dists[i].begin()->first).c_str());
-                               remove((dists[i].begin()->second).c_str());
+                               m->mothurRemove((dists[i].begin()->first));
+                               m->mothurRemove((dists[i].begin()->second));
                        }
                        dists.clear();
                }
@@ -639,7 +749,7 @@ int SplitMatrix::splitDistanceRAM(){
                int numGroups = 0;
 
                ifstream dFile;
-               openInputFile(distFile, dFile);
+               m->openInputFile(distFile, dFile);
 
                while(dFile){
                        string seqA, seqB;
@@ -647,7 +757,7 @@ int SplitMatrix::splitDistanceRAM(){
 
                        dFile >> seqA >> seqB >> dist;
                        
-                       if (m->control_pressed) {   dFile.close();  for(int i=0;i<numGroups;i++){       if(groups[i].size() > 0){  remove((distFile + "." + toString(i) + ".temp").c_str()); }  } return 0; }
+                       if (m->control_pressed) {   dFile.close();  for(int i=0;i<numGroups;i++){       if(groups[i].size() > 0){  m->mothurRemove((distFile + "." + toString(i) + ".temp")); }  } return 0; }
                                        
                        if(dist < cutoff){
                                //cout << "in cutoff: " << dist << endl;
@@ -722,7 +832,7 @@ int SplitMatrix::splitDistanceRAM(){
                                        }
                                }
                        }
-                       gobble(dFile);
+                       m->gobble(dFile);
                }
                dFile.close();