From: Tim Reddy Tim Date: Tue, 23 Dec 2008 19:26:07 +0000 (+0000) Subject: Working on a bug fix... X-Git-Url: http://woldlab.caltech.edu/gitweb/?p=htsworkflow.git;a=commitdiff_plain;h=b1cccd17cfa91f3d023c5469a802f1a21ed0f838 Working on a bug fix... --- diff --git a/htswanalysis/src/GetReadsInSnps/getRegionConsensus.cpp b/htswanalysis/src/GetReadsInSnps/getRegionConsensus.cpp index 1f0ae1b..4cbd878 100755 --- a/htswanalysis/src/GetReadsInSnps/getRegionConsensus.cpp +++ b/htswanalysis/src/GetReadsInSnps/getRegionConsensus.cpp @@ -169,6 +169,8 @@ class Nuc { } char consensus() { + unsigned int N = size(); + if(N == 0) { return ' '; } unsigned int max = 0; unsigned int max_idx = 0; for(unsigned int i = 0; i < 4; i++) { if(n[i] > max) { max = n[i]; max_idx = i; } } @@ -195,8 +197,7 @@ class Nuc { } } - unsigned int N = size(); - if(N == 0) { return ' '; } else if(N < 10) { return tolower(c); } else { return c; } + if(N < 10) { return tolower(c); } else { return c; } } }; @@ -222,6 +223,8 @@ class Window : public Loci { this->name = name; this->length = length; this->sequence = ""; + + if(length > 10000) { cerr << "ERROR: window of size " << length << endl; exit(1); } seq.resize(length); this->reads = 0; @@ -240,6 +243,7 @@ class Window : public Loci { } Window& operator=(const Window& r) { + if (this == &r) return *this; Loci::operator=(r); this->name = r.name; this->length = r.length; @@ -266,6 +270,7 @@ class Window : public Loci { this->reads++; for(unsigned int i = 0; i < r.length(); i++) { int seq_idx = offset + i; + if(seq_idx > 0 && (unsigned)seq_idx >= seq.size()) { cerr << "Error: offset: " << offset << " seq_idx " << seq_idx << " size: " << seq.size() << " seqlen: " << sequence.length() << endl; } if(seq_idx < 0 || (seq_idx >= 0 && (unsigned)seq_idx > this->length) ) { continue; } seq[offset + i].add_nuc(r[i]); } @@ -306,13 +311,17 @@ class Window : public Loci { unsigned int max_len = sequence.length() - offset; unsigned int len = (line_len > max_len)?max_len:line_len; for(unsigned int i = offset; i < offset+len; i++) { + if(i >= seq.size()) { + cerr << "Error: offset of " << i << " exceeds read data: " << seq.size() << " in string of " << sequence.length() << endl; + exit(1); + } char con = toupper(seq[i].consensus()); // weak consensus if lowercase. bool weak_con = seq[i].consensus() != con; if(con == ' ' || weak_con || toupper(con) == toupper(sequence[i])) { - output += sequence[i]; + output.append(1,sequence[i]); } else { - output += con; + output.append(1,con); char buff[128]; sprintf(buff,"%d:%c>%c",i,sequence[i],con); string var = buff; @@ -716,7 +725,7 @@ void count_read_in_features(Windows& windows, Reads& data) { //stop if we have run out of features. if(wind_it == windows.end()) { break; } - if(i->pos + i->length > wind_it->pos && i->pos < (wind_it->pos + wind_it->length)) { + if(i->pos + i->length() > wind_it->pos && i->pos < (wind_it->pos + wind_it->length)) { wind_it->add_read(*i); } } @@ -739,7 +748,7 @@ void retrieveSequenceData(ChromList chrom_filenames, Windows& peaks) { offset = chrom_file.gcount(); } unsigned int begin = i->pos - 1; - unsigned int end = i->pos+i->length-2; + unsigned int end = i->pos+i->length; unsigned int begin_pos = offset + (int)begin/50 + begin; unsigned int end_pos = offset + (int)end/50 + end;