flowdata.cpp

   1 /*
   2  *  flowdata.cpp
   3  *  Mothur
   4  *
   5  *  Created by Pat Schloss on 12/22/10.
   6  *  Copyright 2010 Schloss Lab. All rights reserved.
   7  *
   8  */
   9
  10 #include "flowdata.h"
  11
  12 //**********************************************************************************************************************
  13
  14 FlowData::FlowData(){}
  15
  16 //**********************************************************************************************************************
  17
  18 FlowData::~FlowData(){  /*      do nothing      */      }
  19
  20 //**********************************************************************************************************************
  21
  22 FlowData::FlowData(int numFlows, float signal, float noise, int maxHomoP, string baseFlow) :
  23                         numFlows(numFlows), signalIntensity(signal), noiseIntensity(noise), maxHomoP(maxHomoP), baseFlow(baseFlow){
  24
  25         try {
  26                 m = MothurOut::getInstance();
  27
  28                 flowData.assign(numFlows, 0);
  29 //              baseFlow = "TACG";
  30                 seqName = "";
  31                 locationString = "";
  32         }
  33         catch(exception& e) {
  34                 m->errorOut(e, "FlowData", "FlowData");
  35                 exit(1);
  36         }
  37
  38 }
  39
  40 //**********************************************************************************************************************
  41
  42 bool FlowData::getNext(ifstream& flowFile){
  43
  44         try {
  45                 flowFile >> seqName >> endFlow;
  46         if (seqName.length() != 0) {
  47             //cout << "in Flowdata " + seqName << endl;
  48             for(int i=0;i<numFlows;i++) {       flowFile >> flowData[i];        }
  49             //cout << "in Flowdata read " << seqName + " done" << endl;
  50             updateEndFlow();
  51             translateFlow();
  52             m->gobble(flowFile);
  53                 }else{ m->mothurOut("Error in reading your flowfile, at position " + toString(flowFile.tellg()) + ". Blank name."); m->mothurOutEndLine(); }
  54
  55                 if(flowFile){   return 1;       }
  56                 else            {       return 0;       }
  57         }
  58         catch(exception& e) {
  59                 m->errorOut(e, "FlowData", "getNext");
  60                 exit(1);
  61         }
  62
  63 }
  64
  65 //**********************************************************************************************************************
  66
  67 void FlowData::updateEndFlow(){
  68         try{
  69
  70                 //int currLength = 0;
  71                 float maxIntensity = (float) maxHomoP + 0.49;
  72
  73                 int deadSpot = 0;
  74
  75                 while(deadSpot < endFlow){
  76                         int signal = 0;
  77                         int noise = 0;
  78
  79                         for(int i=0;i<4;i++){
  80                                 float intensity = flowData[i + deadSpot];
  81                                 if(intensity > signalIntensity){
  82                                         signal++;
  83
  84                                         if(intensity  < noiseIntensity || intensity > maxIntensity){
  85                                                 noise++;
  86                                         }
  87                                 }
  88                         }
  89
  90                         if(noise > 0 || signal == 0){
  91                                 break;
  92                         }
  93
  94                         deadSpot += 4;
  95                 }
  96                 endFlow = deadSpot;
  97
  98         }
  99         catch(exception& e) {
 100                 m->errorOut(e, "FlowData", "findDeadSpot");
 101                 exit(1);
 102         }
 103 }
 104
 105 //**********************************************************************************************************************
 106
 107 void FlowData::translateFlow(){
 108
 109         try{
 110                 sequence = "";
 111                 for(int i=0;i<endFlow;i++){
 112                         int intensity = (int)(flowData[i] + 0.5);
 113                         char base = baseFlow[i % 4];
 114
 115                         for(int j=0;j<intensity;j++){
 116                                 sequence += base;
 117                         }
 118                 }
 119
 120                 if(sequence.size() > 4){
 121                         sequence = sequence.substr(4);
 122                 }
 123                 else{
 124                         sequence = "NNNN";
 125                 }
 126         }
 127         catch(exception& e) {
 128                 m->errorOut(e, "FlowData", "translateFlow");
 129                 exit(1);
 130         }
 131 }
 132
 133 //**********************************************************************************************************************
 134
 135 void FlowData::capFlows(int mF){
 136
 137         try{
 138
 139                 maxFlows = mF;
 140                 if(endFlow > maxFlows){ endFlow = maxFlows;     }
 141         translateFlow();
 142
 143         }
 144         catch(exception& e) {
 145                 m->errorOut(e, "FlowData", "capFlows");
 146                 exit(1);
 147         }
 148 }
 149
 150 //**********************************************************************************************************************
 151
 152 bool FlowData::hasMinFlows(int minFlows){
 153
 154         try{
 155                 bool pastMin = 0;
 156                 if(endFlow >= minFlows){        pastMin = 1;    }
 157
 158                 return pastMin;
 159         }
 160         catch(exception& e) {
 161                 m->errorOut(e, "FlowData", "hasMinFlows");
 162                 exit(1);
 163         }
 164 }
 165
 166 //**********************************************************************************************************************
 167
 168 Sequence FlowData::getSequence(){
 169
 170         try{
 171                 return Sequence(seqName, sequence);
 172         }
 173         catch(exception& e) {
 174                 m->errorOut(e, "FlowData", "getSequence");
 175                 exit(1);
 176         }
 177 }
 178
 179 //**********************************************************************************************************************
 180
 181 void FlowData::printFlows(ofstream& outFlowFile){
 182         try{
 183 //      outFlowFile << '>' << seqName << locationString << " length=" << seqLength << " numflows=" << maxFlows << endl;
 184                 outFlowFile << seqName << ' ' << endFlow << ' ' << setprecision(2);
 185
 186                 for(int i=0;i<maxFlows;i++){
 187                         outFlowFile << flowData[i] << ' ';
 188                 }
 189                 outFlowFile << endl;
 190         }
 191         catch(exception& e) {
 192                 m->errorOut(e, "FlowData", "printFlows");
 193                 exit(1);
 194         }
 195 }
 196
 197 //**********************************************************************************************************************
 198
 199 void FlowData::printFlows(ofstream& outFlowFile, string scrapCode){
 200         try{
 201                 outFlowFile << seqName << '|' << scrapCode << ' ' << endFlow << ' ' << setprecision(2);
 202
 203                 for(int i=0;i<numFlows;i++){
 204                         outFlowFile << flowData[i] << ' ';
 205                 }
 206                 outFlowFile << endl;
 207         }
 208         catch(exception& e) {
 209                 m->errorOut(e, "FlowData", "printFlows");
 210                 exit(1);
 211         }
 212 }
 213
 214 //**********************************************************************************************************************
 215
 216 string FlowData::getName(){
 217
 218         try{
 219                 return seqName;
 220         }
 221         catch(exception& e) {
 222                 m->errorOut(e, "FlowData", "getName");
 223                 exit(1);
 224         }
 225 }
 226
 227 //**********************************************************************************************************************