added shhh.seqs command

[mothur.git] / sffinfocommand.cpp

diff --git a/sffinfocommand.cpp b/sffinfocommand.cpp
index d0ec9e2a93b5109c91e4a940f8b3bbaa64bdfe41..3562dcc3e88265a93c23af38fa69b3b2c776e4c1 100644 (file)
--- a/sffinfocommand.cpp
+++ b/sffinfocommand.cpp
@@ -316,7 +316,7 @@ int SffInfoCommand::extractSffInfo(string input, string accnos){
                if (sfftxt) { m->openOutputFile(sfftxtFileName, outSfftxt); outSfftxt.setf(ios::fixed, ios::floatfield); outSfftxt.setf(ios::showpoint);  outputNames.push_back(sfftxtFileName);  outputTypes["sfftxt"].push_back(sfftxtFileName); }
                if (fasta)      { m->openOutputFile(outFastaFileName, outFasta);        outputNames.push_back(outFastaFileName); outputTypes["fasta"].push_back(outFastaFileName); }
                if (qual)       { m->openOutputFile(outQualFileName, outQual);          outputNames.push_back(outQualFileName); outputTypes["qual"].push_back(outQualFileName);  }
-               if (flow)       { m->openOutputFile(outFlowFileName, outFlow);          outputNames.push_back(outFlowFileName);  outputTypes["flow"].push_back(outFlowFileName);  }
+               if (flow)       { m->openOutputFile(outFlowFileName, outFlow);          outputNames.push_back(outFlowFileName);  outFlow.setf(ios::fixed, ios::floatfield); outFlow.setf(ios::showpoint); outputTypes["flow"].push_back(outFlowFileName);  }
                
                ifstream in;
                in.open(input.c_str(), ios::binary);
@@ -331,8 +331,9 @@ int SffInfoCommand::extractSffInfo(string input, string accnos){
                if (header.version != "0001") { m->mothurOut("Version is not supported, only support version 0001."); m->mothurOutEndLine(); return count; }
        
                //print common header
-               if (sfftxt) { printCommonHeader(outSfftxt, header); }
-       
+               if (sfftxt) {   printCommonHeader(outSfftxt, header);           }
+               if (flow)       {       outFlow << header.numFlowsPerRead << endl;      }
+                       
                //read through the sff file
                while (!in.eof()) {
                        
@@ -515,6 +516,9 @@ int SffInfoCommand::readHeader(ifstream& in, Header& header){
                        if (header.name.length() > header.nameLength) { header.name = header.name.substr(0, header.nameLength);  }
                        delete[] tempBuffer;
                        
+                       //extract info from name
+                       decodeName(header.timestamp, header.region, header.xy, header.name);
+                       
                        /* Pad to 8 chars */
                        unsigned long int spotInFile = in.tellg();
                        unsigned long int spot = (spotInFile + 7)& ~7;
@@ -585,6 +589,43 @@ int SffInfoCommand::readSeqData(ifstream& in, seqRead& read, int numFlowReads, i
        }
 }
 //**********************************************************************************************************************
+int SffInfoCommand::decodeName(string& timestamp, string& region, string& xy, string name) {
+       try {
+               
+               string time = name.substr(0, 6);
+               unsigned int timeNum = m->fromBase36(time);
+                       
+               int q1 = timeNum / 60;
+               int sec = timeNum - 60 * q1;
+               int q2 = q1 / 60;
+               int minute = q1 - 60 * q2;
+               int q3 = q2 / 24;
+               int hr = q2 - 24 * q3;
+               int q4 = q3 / 32;
+               int day = q3 - 32 * q4;
+               int q5 = q4 / 13;
+               int mon = q4 - 13 * q5;
+               int year = 2000 + q5;
+               
+               timestamp = toString(year) + "_" + toString(mon) + "_" + toString(day) + "_" + toString(hr) + "_" + toString(minute) + "_" + toString(sec);
+               
+               region = name.substr(7, 2);
+               
+               string xyNum = name.substr(9);
+               unsigned int myXy = m->fromBase36(xyNum);
+               int x = myXy >> 12;
+               int y = myXy & 4095;
+               
+               xy = toString(x) + "_" + toString(y);
+                       
+               return 0;
+       }
+       catch(exception& e) {
+               m->errorOut(e, "SffInfoCommand", "decodeName");
+               exit(1);
+       }
+}
+//**********************************************************************************************************************
 int SffInfoCommand::printCommonHeader(ofstream& out, CommonHeader& header) {
        try {
        
@@ -612,9 +653,9 @@ int SffInfoCommand::printHeader(ofstream& out, Header& header) {
        try {
                
                out << ">" << header.name << endl;
-               out << "Run Prefix: " << endl;
-               out << "Region #:  " << endl;
-               out << "XY Location: " << endl << endl;
+               out << "Run Prefix: " << header.timestamp << endl;
+               out << "Region #:  " << header.region << endl;
+               out << "XY Location: " << header.xy << endl << endl;
                
                out << "Run Name:  " << endl;
                out << "Analysis Name:  " << endl;
@@ -693,7 +734,7 @@ int SffInfoCommand::printFastaSeqData(ofstream& out, seqRead& read, Header& head
                        }
                }
                
-               out << ">" << header.name << endl;
+               out << ">" << header.name  << " xy=" << header.xy << endl;
                out << seq << endl;
                
                return 0;
@@ -713,15 +754,15 @@ int SffInfoCommand::printQualSeqData(ofstream& out, seqRead& read, Header& heade
                                out << "0\t0\t0\t0";
                        }
                        else if((header.clipQualRight != 0) && ((header.clipQualRight-header.clipQualLeft) >= 0)){
-                               out << ">" << header.name << " length=" << (header.clipQualRight-header.clipQualLeft) << endl;
+                               out << ">" << header.name << " xy=" << header.xy << " length=" << (header.clipQualRight-header.clipQualLeft) << endl;
                                for (int i = (header.clipQualLeft-1); i < (header.clipQualRight-1); i++) {   out << read.qualScores[i] << '\t'; }
                        }
                        else{
-                               out << ">" << header.name << " length=" << (header.clipQualRight-header.clipQualLeft) << endl;
+                               out << ">" << header.name << " xy=" << header.xy << " length=" << (header.clipQualRight-header.clipQualLeft) << endl;
                                for (int i = (header.clipQualLeft-1); i < read.qualScores.size(); i++) {   out << read.qualScores[i] << '\t';   }                       
                        }
                }else{
-                       out << ">" << header.name << " length=" << read.qualScores.size() << endl;
+                       out << ">" << header.name << " xy=" << header.xy << " length=" << read.qualScores.size() << endl;
                        for (int i = 0; i < read.qualScores.size(); i++) {   out << read.qualScores[i] << '\t';  }
                }
                
@@ -738,10 +779,16 @@ int SffInfoCommand::printQualSeqData(ofstream& out, seqRead& read, Header& heade
 //**********************************************************************************************************************
 int SffInfoCommand::printFlowSeqData(ofstream& out, seqRead& read, Header& header) {
        try {
+               if(header.clipQualRight > header.clipQualLeft){
+                       
+                       int rightIndex = 0;
+                       for (int i = 0; i < header.clipQualRight; i++) {  rightIndex +=  read.flowIndex[i];     }
+
+                       out << header.name << ' ' << rightIndex;
+                       for (int i = 0; i < read.flowgram.size(); i++) { out << setprecision(2) << ' ' << (read.flowgram[i]/(float)100);  }
+                       out << endl;
+               }
                
-               out << ">" << header.name << endl;
-               for (int i = 0; i < read.flowgram.size(); i++) { out << setprecision(2) << (read.flowgram[i]/(float)100) << '\t';  }
-               out << endl;
                
                return 0;
        }