Skip to content

Commit

Permalink
src/load_data_for_complexity.cpp: fixing a bug in checking whether a …
Browse files Browse the repository at this point in the history
…read is mapped -- the condition was inverted
  • Loading branch information
andrewdavidsmith committed Oct 19, 2024
1 parent 0a94079 commit 251d39f
Showing 1 changed file with 11 additions and 6 deletions.
17 changes: 11 additions & 6 deletions src/load_data_for_complexity.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -520,6 +520,11 @@ load_coverage_counts_GR(const string &input_file_name, const uint64_t seed,
#ifdef HAVE_HTSLIB
// Deal with SAM/BAM format only if we have htslib

static inline void
is_unmapped(const bamxx::bam_rec &aln) {
return get_tid(aln) == -1;
}

static inline void
swap(bamxx::bam_rec &a, bamxx::bam_rec &b) {
std::swap(a.b, b.b);
Expand Down Expand Up @@ -593,12 +598,12 @@ load_counts_BAM(const uint32_t n_threads, const string &inputfile,

// find first mapped read to start
bamxx::bam_rec aln;
while (hts.read(hdr, aln) && get_tid(aln) == -1)
while (hts.read(hdr, aln) && is_unmapped(aln))
;

size_t n_reads{};
// if all reads unmapped, must return
if (get_tid(aln) == -1)
if (is_unmapped(aln))
return n_reads;

// to check that reads are sorted properly
Expand All @@ -611,7 +616,7 @@ load_counts_BAM(const uint32_t n_threads, const string &inputfile,
size_t current_count = 1;

while (hts.read(hdr, aln)) {
if (get_tid(aln) != -1)
if (is_unmapped(aln))
continue; // skip unmapped reads

const aln_pos_t curr{aln};
Expand Down Expand Up @@ -739,11 +744,11 @@ load_coverage_counts_BAM(const uint32_t n_threads, const string &inputfile,

// find first mapped read to start
bamxx::bam_rec aln;
while (hts.read(hdr, aln) && get_tid(aln) == -1)
while (hts.read(hdr, aln) && is_unmapped(aln))
;

size_t n_reads{};
if (get_tid(aln) == -1) // no reads unmapped
if (is_unmapped(aln)) // no reads unmapped
return 0;

// to check reads are sorted properly
Expand All @@ -761,7 +766,7 @@ load_coverage_counts_BAM(const uint32_t n_threads, const string &inputfile,
const hts_pos_t max_dist = bin_size + max_width;

while (hts.read(hdr, aln)) {
if (get_tid(aln) == -1)
if (is_unmapped(aln))
continue; // check that read is mapped

const hts_pos_t rlen = rlen_from_cigar(aln);
Expand Down

0 comments on commit 251d39f

Please sign in to comment.